RF-Track Reference Manual

RF-Track Reference Manual Version 2.5.4 Andrea Latina AUTHOR AND CONTACT: Andrea Latina Beams Department Accelerator and Beam Physics Group CERN CH-1211 GENEVA 23 SWITZERLAND [email protected] Copyright © 2025 by CERN, Geneva, Switzerland Copyright and any other appropriate legal protection of this computer program and associated documentation reserved in all countries of the world. Organisations collaborating with CERN may receive this program and documentation freely and without charge. CERN undertakes no obligation for the maintenance of this program, nor responsibility for its correctness, and accepts no liability whatsoever resulting from its use. Program and documentation are provided solely for the use of the organisation to which they are distributed. This program may not be copied or otherwise distributed without permission. This message must be retained on this and any other authorised copies. The material cannot be sold. CERN should be given credit in all references. Started in April 2020 Contents IUser Manual 1Introduction .................................................. 11 1.1 Getting started 11 1.1.1 Conventions used in this manual . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 12 1.2 Installing RF-Track 12 1.2.1 OnmacOS..................................................... 12 1.2.2 OnLinux....................................................... 12 1.2.3 OnWindows.................................................... 13 1.3 Running RF-Track 13 1.3.1 Anexampleprogram ............................................ 13 1.3.2 Physicalunits ................................................... 16 1.3.3 Predefinedconstants............................................. 17 1.4 Run-time parameters and options 17 1.4.1 Parallelism ..................................................... 17 1.4.2 Environmentvariables ............................................ 18 1.4.3 RandomNumberGenerator ....................................... 18 1.4.4 Versionnumber ................................................. 18 1.5 Notes for Python users 19 1.6 Further information 19 1.6.1 Acknowledgments............................................... 19 1.6.2 CitingRF-Trackinpublications ...................................... 20 2Beam Models ................................................ 21 2.1 Tracking in time or in space 21 2.2 Single-bunch beams 21 2.2.1 Bunch6d....................................................... 21 2.2.2 Bunch6dT...................................................... 25 2.2.3 Conversion between Bunch6d and Bunch6dT . . . . . . . . . . . . . . . . . . . . . . . . . 28 2.2.4 Coastingbeams ................................................ 28 2.2.5 Particleslifetime................................................. 28 2.2.6 Offsettingthebeam ............................................. 29 2.3 Multi-bunch beams 29 2.3.1 ExampleofBeamDefinition ....................................... 29 2.4 Twiss parameters 31 2.4.1 Bunch6d_twiss .................................................. 31 2.4.2 Bunch6dT_twiss ................................................. 31 2.5 Inquiring the bunch properties 32 2.5.1 Bunch6d_info................................................... 32 2.5.2 Bunch6dT_info .................................................. 33 2.6 Bunch persistency 35 2.6.1 Saving and loading a beam on disk . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 35 2.6.2 ExportingasDSTfile .............................................. 35 2.6.3 ExportingasSDDSfile............................................. 35 2.6.4 Savingthephasespace .......................................... 35 2.7 Spin Polarization Tracking 36 2.7.1 Settingthepolarization ........................................... 36 2.7.2 Inquiringthepolarization.......................................... 37 2.7.3 Example of spin polarization tracking . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 37 3Beam Tracking ............................................... 41 3.1 Tracking environments 41 3.2 Lattice 42 3.2.1 Addingelements ................................................ 42 3.2.2 ImportingaMAD-Xlattice......................................... 42 3.2.3 Elementmisalignment ............................................ 43 3.2.4 Trackingthebeam............................................... 44 3.2.5 Accessing and modifying the elements . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 44 3.3 Volume 45 3.3.1 Addingelements ................................................ 45 3.3.2 Accessing and modifying the elements . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 46 3.3.3 VolumeasaLatticeelement....................................... 51 3.3.4 TransporttableandScreens ....................................... 51 3.4 Particle losses 52 3.5 Synchronization of time-dependent elements with the beam 52 3.5.1 Autophasing ................................................... 53 3.5.2 Automatic setting of magnetic strengths . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 54 3.6 Backtracking 54 4Beamline Elements ........................................... 55 4.1 Introduction 55 4.1.1 Methods available to all elements . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 55 4.1.2 Inquiring and plotting the electromagnetic field . . . . . . . . . . . . . . . . . . . . . . . 56 4.1.3 Trackingthebeamproperties ...................................... 56 4.2 Matrix-based elements 60 4.2.1 Drift .......................................................... 60 4.2.2 Quadrupole.................................................... 61 4.2.3 Sectorbendingdipole............................................ 63 4.2.4 Rectangularbendingdipole....................................... 65 4.2.5 Corrector ...................................................... 66 4.3 Special elements 67 4.3.1 Coil........................................................... 67 4.3.2 Solenoid....................................................... 68 4.3.3 Undulator...................................................... 69 4.3.4 TransferLine .................................................... 70 4.3.5 Travelling-wavestructure.......................................... 72 4.3.6 Standing-wavestructure .......................................... 74 4.3.7 Pillboxcavity ................................................... 78 4.3.8 Multipolemagnet ............................................... 79 4.3.9 Absorber ...................................................... 81 4.3.10 Electroncooler ................................................. 82 4.3.11 Adiabaticmatchingdevice ....................................... 83 4.3.12 Space-chargeField.............................................. 84 4.3.13 Travelling-WaveField ............................................. 85 4.3.14 LaserBeam..................................................... 86 4.3.15 VolumeasaLatticeelement....................................... 87 4.4 Field maps 88 4.4.1 1DRFfieldmaps ................................................ 89 4.4.2 2DRFfieldmaps ................................................ 91 4.4.3 3DRFfieldmaps ................................................ 93 4.4.4 1Dstaticmagneticfieldmaps...................................... 95 4.4.5 2Dstaticmagneticfieldmaps...................................... 96 4.4.6 3Dstaticmagneticfieldmaps...................................... 97 4.5 Beam diagnostics 98 4.5.1 Beampositionmonitor............................................ 98 4.5.2 Screens ....................................................... 99 5Collective Effects ............................................ 101 5.1 Space charge 103 5.1.1 Space-chargeinVolume......................................... 103 5.1.2 Space-chargeinLattice ......................................... 103 5.1.3 SpaceChargeModels .......................................... 104 5.2 Intra-Beam Scattering 106 5.3 Incoherent synchrotron radiation 107 5.4 Magnetic multipolar errors 107 5.5 Beam loading 108 5.5.1 Beam loading in ultrarelativistic scenarios . . . . . . . . . . . . . . . . . . . . . . . . . . . . 108 5.5.2 Beam loading in standing-wave structures . . . . . . . . . . . . . . . . . . . . . . . . . . . 110 5.6 Wakefields 112 5.6.1 Short-rangewakefield ........................................... 112 5.6.2 Long-rangewakefield ........................................... 114 5.6.3 Genericwakefield .............................................. 116 5.7 Passage of particles through matter 117 5.7.1 MultipleCoulombscattering...................................... 117 5.7.2 Stoppingpower ................................................ 117 5.7.3 Energystraggling ............................................... 118 6Inverse Compton Scattering .................................. 119 6.1 Introduction 119 6.1.1 Generalparameters ............................................ 119 6.1.2 Particles-laserinteraction......................................... 120 6.2 Collecting the generated photons 121 7Bunch Generation at a Photocathode ........................ 123 7.1 Introduction 123 7.2 Generator options 123 7.2.1 Particleproperties .............................................. 123 7.2.2 Longitudinaldistribution.......................................... 123 7.2.3 Transverse spatial distribution . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 124 7.2.4 Transverse momentum distribution . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 124 7.2.5 Additionalparameters........................................... 124 7.3 Example of Generator 124 7.4 Photo emission simulation 125 7.5 Mirror charges at cathode 126 7.5.1 Misalignedcathode ............................................ 126 8Extending RF-Track .......................................... 127 8.1 Introduction 127 8.2 Custom Lattice Elements with UserElement 127 8.2.1 Keyusecases ................................................. 127 8.2.2 Implementation................................................ 128 8.2.3 Example: A Custom Element in Python . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 128 8.3 Custom Electromagnetic Fields with UserField 129 8.3.1 Applications................................................... 129 8.3.2 FieldEvaluation ................................................ 129 8.3.3 Example: A Custom Electromagnetic field in Python . . . . . . . . . . . . . . . . . . . 130 8.4 Custom collective effects with UserEffect 131 8.4.1 Capabilities ................................................... 131 8.4.2 Returningmodifiedbunches ...................................... 131 8.5 Secondary particle generation and external interfaces 131 8.5.1 Example: A Custom Collective Effect in Python . . . . . . . . . . . . . . . . . . . . . . . 131 8.6 Traversing Beamlines with UserVisitor 132 8.6.1 Keyusecases ................................................. 132 8.6.2 Implementing a Visitor in Python . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 133 8.7 Getting started 134 9Examples ................................................... 135 9.1 Example of bunch creation 135 9.1.1 Bunch6d from an arbitrary user-defined distribution . . . . . . . . . . . . . . . . . . . . 135 9.1.2 Chirped Bunch6d from Twiss parameters . . . . . . . . . . . . . . . . . . . . . . . . . . . . 136 Index ....................................................... 137 I 1Introduction ......................... 11 2Beam Models ....................... 21 3Beam Tracking ...................... 41 4Beamline Elements .................. 55 5Collective Effects ................... 101 6Inverse Compton Scattering ......... 119 7Bunch Generation at a Photocathode 123 8Extending RF-Track ................. 127 9Examples .......................... 135 Index .............................. 137 User Manual 16 Chapter 1. Introduction Line 47 tracks the bunch B0 through the FODO cell and stores the outcoming bunch in B1 , another instance of type Bunch6d(). Lines 50-51 retrieve the Twiss parameters and the phase space and store them in two Octave variables Tand M. Lines 54-65 plot the results using the standard Octave plotting commands. Figure 1.1 shows the output of the example: the Twiss parameters, βx and βy , and the phase space plot. Figure 1.1: The output of the example. 1.3.2 Physical units Internally, RF-Track stores each physical quantity in the most suitable units for numerical computation in accelerator physics. Table 1.1 shows the units grouped by conceptual category. Regarding positions and angles, notice that all beam-related quantities are in units of millimeters or milliradians, whereas all machine-related quantities are expressed in meters and radians. Table 1.1: RF-Track physical units. Quantity Symbols Unit Bunch population Nnumber of particles Particle mass mMeV/c2 Particle charge Q e Particle positions x,y,zmm Particle angles x0,y0mrad Particle momenta Px,Py,Pz,PMeV/c Particle energy EMeV Time tmm/c Element offsets and positions Xo,Yo,Zom Element pitch X0 orad Element yaw Y0 orad Element roll Z0 orad 1.4 Run-time parameters and options 17 1.3.3 Predefined constants Several constants are predefined within RF-Track for the user’s convenience: clight % speed of light, in m/s muonmass % muon mass in MeV/c^2 protonmass % proton mass, in MeV/c^2 electronmass % electron mass, in MeV/c^2 muonlifetime % muon lifetime, in mm/c s, ms, us, ns, ps, fs % various units if time, in mm/c C, mC, uC, nC, pC % various units of charge, in e For example, if one needs to define a time interval of 5 ps, say dt, one can write: octave:1> dt = 5 *RF_Track.ps dt = 1.4990 % 5 ps in mm/c 1.4 Run-time parameters and options Once RF-Track is loaded, a few variables become available to the user to customize how RF-Track operates. 1.4.1 Parallelism RF-Track is a parallel application. By default, RF-Track uses the maximum number of parallel threads available on your machine, which can be inquired via the RF-Track’s variable “ max_number_of_threads ”. The user can change this number (typically reducing it to avoid overloading the CPU) by setting the variable number_of_threads. In Octave, % The number of threads that RF-Track must use RF_Track.number_of_threads = <YOUCHOOSE>; % The max number of threads supported by your CPU [read-only] RF_Track.max_number_of_threads In Python, % The number of threads that RF-Track must use RF_Track.cvar.number_of_threads = <YOUCHOOSE> % The max number of threads supported by your CPU [read-only] RF_Track.max_number_of_threads The desired number of threads that RF-Track should use (in this example as <YOUCHOOSE> ) should not exceed RF_Track.max_number_of_threads . Note that the fastest RF-Track setup isn’t necessarily achieved by using all the threads available on your system. Sometimes, using less can 18 Chapter 1. Introduction lead to faster simulations. This depends on the CPU you are using, the number of particles you are tracking, and the simulation setup in terms of collective effects, transport table, and screens, i.e., the points where RF-Track needs to break parallel tracking to compute average or global bunch quantities. Each simulation case should be studied individually and a preliminary test with different numbers of threads is recommended to find the optimum number of threads. 1.4.2 Environment variables A number of environment variables exist to customize RF-Track’s runtime behaviour: RF_TRACK_NUMBER_OF_THREADS The number of threads that RF-Track uses can also be specified via the environment variable, RF_TRACK_NUMBER_OF_THREADS . This allows the user to specify the number of threads at runtime without modifying the input script. In this example: $ RF_TRACK_NUMBER_OF_THREADS=8 octave-cli my_rft_script.m Octave will execute the RF-Track script my_rft_script.m using a maximum of eight parallel threads. Similarly, with Python: $ RF_TRACK_NUMBER_OF_THREADS=8 python my_rft_script.py RF_TRACK_NO_UPDATE_CHECK At start-up, RF-Track checks if a newer version is available. To turn off the automatic check, set this environment variable to ’1’. Example: $ export RF_TRACK_NO_UPDATE_CHECK=1 $ python my_rft_script.py 1.4.3 Random Number Generator RF-Track uses high-quality random number generators. The user can use the following commands to set the desired random number generator and its initial seed. rng_set(name) % change the random number generator (RNG) rng_set_seed(seed) % set the seed for the RNG rng_get() % return what is the current RNG The default random number generator is “mt19937” of Makoto Matsumoto and Takuji Nishimura, known as the “Mersenne Twister” generator, which is among the fastest high-quality generators available. Table 1.2 shows a list of the random number generators available to RF-Track. Consult [3] for a detailed description of each of them. 1.4.4 Version number RF_Track.version % the RF-Track version number [read-only] 1.5 Notes for Python users 19 Name Description taus2 Maximally equidistributed combined Tausworthe generator by L’Ecuyer mt19937 Makoto Matsumoto and Takuji Nishimura generator gfsr4 Lagged-fibonacci generator ranlxs0 Second-generation version of the RANLUX algorithm of Luscher ranlxs1 Like the previous by with increased order of strength ranlxs2 Like the previous by with increased order of strength mrg Fifth-order multiple-recursive generator by L’Ecuyer, Blouin and Coutre ranlux Implementation of the original algorithm developed by Luscher ranlux389 Like the previous but gives the highest level of randomness ranlxd1 Double precision output (48 bits) from the RANLXS generator ranlxd2 Like the previous by with increased order of strength Table 1.2: List of random number generators available to RF-Track Returns the RF-Track version number. 1.5 Notes for Python users Three important points for Python users: 1. In RF-Track, the input arguments of each function are positional, not keyword arguments as Python users may be used to and expect. This means that you cannot change their order. 2. In this manual, the syntax function(arg1=var1, arg2=var2) is purely illustrative, as it is intended to help the reader immediately know the meaning of each input parameter and the default value if the parameter is omitted. It does not indicate the use of keyword arguments. 3. It is also important to know that whenever RF-Track’s inputs are vectors or matrices, RFTrack expects numpy arrays. Therefore, Octave lines like Quad = Multipole (length, [0, k1L]) Sext = Multipole (length, [0, 0, k2L]) must be written in Python as import RF_Track as rft import numpy as np Quad = rft.Multipole(length, np.array([0, k1L])) Sext = rft.Multipole(length, np.array([0, 0, k2L])) 1.6 Further information 1.6.1 Acknowledgments The author gratefully acknowledges the many colleagues who contributed to the development and validation of RF-Track. Early users — Stefano Benedetti, Costanza Agazzi, Yanliang Han, Yongke 20 Chapter 1. Introduction Zhao, Luca Garolfi, Elena Fol, Alexander Malyzhenkov, and Vlad Musat — provided essential feedback that strengthened the code and guided its adaptation to real-life applications. Feature contributions came from Mohsen Dayyani Kelisani (analytic RF structures), Bernd Michael Stechauner (multiple Coulomb scattering models), Avni Aksoy (ASTRA-like particle generator), and Riccardo De Maria (TPSA clarifications). PhD students Javier Olivares Herrador and Paula Desiré Valdor contributed the beam-loading and intra-beam scattering models, while Summer Student Lale Balat helped build a robust test suite. Stewart Boogert interfaced BDSIM to RF-Track. The author also acknowledges collaborations with the CLEAR facility at CERN, the Oncology Department of Lausanne University Hospital (CHUV), and the CLIC and FCC study groups, which provided valuable opportunities for benchmarking and applications, as well as the support of the Accelerator and Beam Physics Group, Beams Department, CERN. 1.6.2 Citing RF-Track in publications We have invested a lot of time and effort in creating RF-Track; please cite it when using it. To cite RF-Track in publications, use: Andrea Latina "RF-Track Reference Manual", CERN, Geneva, Switzerland, 2024 DOI: 10.5281/zenodo.3887085 A BibTeX entry for LaTeX users is: @techreport{, address = {Geneva, Switzerland}, author = {Latina, Andrea}, doi = {10.5281/zenodo.3887085}, institution = {CERN}, title = {RF-Track Reference Manual}, year = {2024} } 2. Beam Models 2.1 Tracking in time or in space RF-Track implements two particle tracking methods: tracking in time and tracking in space. The tracking in time should be preferred in space-charge-dominated regimes, where the relative positions of the particles in space matter. The tracking in space suits better space-charge-free regions, where particles are independent of each other, and they can be transported simultaneously from the entrance plane of an element to its end, element by element. RF-Track provides two distinct beam types to implement these two models: Bunch6dT for tracking in time and Bunch6d for tracking in space. A dedicated tracking environment exists for each of these two beam types: Lattice for Bunch6d , and Volume Bunch6dT . This chapter will describe them in detail. 2.2 Single-bunch beams 2.2.1 Bunch6d When Bunch6d is used, the tracking is performed using the accelerator longitudinal coordinate S as the integration variable. This corresponds to what is described in accelerator physics textbooks, which is at the basis of matrix-based beam optics. In this model, all particles are on the same plane at a given longitudinal coordinate S at the beginning of an element and are then transported to the end of the element, updating the arrival time as the longitudinal coordinate. In Bunch6d , the beam is represented by a set of macro-particles whose state vector is an extended trace space: x,x0,y,y0,t,P,m,Q,N The meaning of each symbol, along with the units used by RF-Track internally, follows: 22 Chapter 2. Beam Models x,ytransverse coordinates [mm] x0,y0transverse angles [mrad] tarrival time at S[mm/c] Ptotal momentum [MeV/c] mmass [MeV/c2] Qcharge of the single particle [e] Nnumber of single particles in each macro-particle [#] Bunch6d also stores the longitudinal coordinate S of the bunch along the Lattice, which can be accessed like this: current_S = B.S; % m  Emittance calculation When using a Bunch6d, the normalised emittances are computed as: εx=βavg γavg detcovx,x01/2 εy=βavg γavg detcovy,y01/2 εz=1 mdet(cov{t,E})1/2 ε4D =βavg γavg detcovx,x0,y,y01/4 ε6D =βavg γavg detcovx,x0,y,y0,t,E Pavg 1/6 Constructors A new Bunch6d can be created in multiple ways. The most common ones are two: from a set of Twiss parameters or directly from a beam matrix with the phase space. Multi-species bunches can be created using a beam matrix of the extended phase space. Here is the list of constructors: B = Bunch6d(mass, population, charge, Pref, Twiss, nParticles, sigmaCut=0 ); B = Bunch6d(mass, population, charge, [ X XP Y YP T P ID ] ); B = Bunch6d(mass, population, charge, [ X XP Y YP T P ] ); B = Bunch6d( [ X XP Y YP T P MASS Q N ID ] ); B = Bunch6d( [ X XP Y YP T P MASS Q N ] ); B = Bunch6d( [ X XP Y YP T P MASS Q ] ); B = Bunch6d(); The possible arguments are: 2.2 Single-bunch beams 23 mass the mass of the single particle [MeV/c2] population the total number of real particles in the bunch [#] charge charge of the single particle [e] Pref reference momentum [MeV/c] Twiss instance of object Bunch6d_twiss (see section 2.4.1) nParticles number of macro-particles in bunch [#] sigmaCut if >0 cuts the distributions at sigma_cut sigmas [#] Xcolumn vector of the horizontal coordinates [mm] Ycolumn vector of the vertical coordinates [mm] Tcolumn vector of the arrival times [mm/c] Pcolumn vector of the total momenta [MeV/c] XP column vector of the x0angles [mrad] YP column vector of the y0angles [mrad] MASS column vector of masses [MeV/c2] Qcolumn vector of single-particle charges [e] N column vector of numbers of single particles per macro particle [#] ID column vector of particle ID’s [INTEGER] With no arguments, the default constructor allows the creation of an empty beam. The return value is an object of type Bunch6d. Accessing the phase space The particle coordinates can be inquired using the method get_phase_space() . This method allows retrieving the bunch’s particle distribution in multiple ways. M = B.get_phase_space(format=’%x %xp %y %yp %t %Pc’,which=’good’); The two arguments are: format_fmt This parameter allows the user to specify and format what phase space representation should be returned which ’all’ | ’good’ . This parameter specified whether all particles should be returned, including the lost ones, or just the good ones (default=’good’) Table 2.1 shows all possible identifiers accepted by get_phase_space() . The return value of this function, M, is a matrix with the requested phase space. Modifying the phase space The particle coordinates can be changed by the user using the method set_phase_space() : B.set_phase_space( [ X XP Y YP T P ] ); The only argument accepted by this method is a 6-column matrix containing the phase space in Bunch6d format, that is, the phase space columns are x , y , the particle’s transverse position in mm, x0 , y0 the angles in mrad, t the arrival time in mm/ c , and P the total momentum in MeV/ c , in the 24 Chapter 2. Beam Models %x horizontal coordinate [mm] %y vertical coordinate [mm] %t arrival time, t[mm/c] %dt relative arrival time, t−t0[mm/c] %z longitudinal coordinate w.r.t. the reference particle [mm] %deg@MHz longitudinal coordinate in degrees at a specified frequency in MHz, e.g. %deg@750 for degrees at 750 MHz [deg] %K kinetic energy [MeV] %E total energy [MeV] %P total momentum [MeV/c] %d relative momentum, δ= (P−P0)/P0[permille] %xp horizontal angle, Px/Pz[mrad] %yp vertical angle, Py/Pz[mrad] %tp inverse of Vz[mrad/c] %Px horizontal momentum, Px[MeV/c] %Py vertical momentum, Py[MeV/c] %Pz longitudinal momentum, Pz[MeV/c] %px normalized horizontal momentum, Px/P0[mrad] %py normalized vertical momentum, Py/P0[mrad] %pz normalized longitudinal momentum, Pz/P0[mrad] %pt normalized relative energy difference, pt= (E−E0)/P0c[permille] %Vx horizontal velocity [c] %Vy vertical velocity [c] %Vz longitudinal velocity [c] %m mass [MeV/c2] %Q charge [e+] %N number of particles per macro particle [#] %id particle’s id [#] Table 2.1: List of identifiers accepted by Bunch6d::get_phase_space(). 2.2 Single-bunch beams 25 order given above (x,x0,y,y0,t,P). 2.2.2 Bunch6dT Bunch6dT allows tracking “in time”. When Bunch6dT is used, the time, t , is the integration variable. The 6D phase-space coordinates of each particle are (X,Y,Z,Px,Py,Pz) . Bunch6dT maintains, t, the clock common to all particles, updated at each integration step. In Bunch6dT , the beam is represented by a set of macro-particles whose state vector is an extended phase space: X,Px,Y,Py,Z,Pz,m,Q,N,t0 The positions and momenta of each particle are updated as they are transported through the accelerator. The meaning of each symbol, along with the units used to store the information internally, follows: X,Y,Ztransverse and longitudinal coordinates [mm] Px,Py,Pztransverse and longitudinal momenta [MeV/c] mmass [MeV/c2] Qcharge of the single particle [e] Nnumber of single particles in each macro-particle [#] t0creation time [mm/c] In Bunch6dT , a major difference from Bunch6d is the presence of t0 , the creation time of each particle. This allows, for example, the simulation of cathodes and particle emission. Bunch6dT also stores the time t at which the bunch is taken, which can be accessed (read and write) like this: current_time = B.t; % mm/c  Emittance calculation When using a Bunch6dT, the normalised emittances are computed as: εx=1 mdet(cov{x,Px})1 2 εy=1 mdet(cov{y,Py})1 2 εz=1 mdet(cov{z,Pz})1 2 ε4D =1 mdet(cov{x,Px,y,Py})1 4 ε6D =1 mdet(cov{x,Px,y,Py,z,Pz})1 6 Constructors A new Bunch6dT can be created in multiple ways. The most common ones are two: from a set of Twiss parameters or directly from a beam matrix with the phase space. Multi-specie bunches can be created using a beam matrix of the extended phase space. Follows the list of constructors: 32 Chapter 2. Beam Models T.disp_py; % rad, vertical dispersion prime T.disp_z; % m, longitudinal dispersion Similarly to Bunch6d_Twiss , the longitudinal phase space can be specified in three alternative ways: 1. Giving both the normalised longitudinal emittance emitt_zand the Twiss parameter βz. 2. Giving the normalised longitudinal emittance emitt_zand either sigma_zor sigma_pz. 3. Giving both sigma_zand sigma_pz. 2.5 Inquiring the bunch properties The user can inquire about a bunch’s statistical quantities, such as the average value or the standard deviation of the phase space variables, the emittances, or the Twiss parameters, by calling the method: I = B.get_info(); Two versions exist, for Bunch6d and for Bunch6dT. 2.5.1 Bunch6d_info If Bis a bunch of type Bunch6d,B.get_info() returns: I = B.get_info(); I.S; % m I.mean_x; % mm, average H position I.mean_y; % mm, average V position I.mean_t; % mm/c, average arrival time I.mean_xp; % mrad, average H angle I.mean_yp; % mrad, average V angle I.mean_Px; % average Px in MeV/c I.mean_Py; % average Py in MeV/c I.mean_Pz; % average Pz in MeV/c I.mean_P; % average momentum in MeV/c I.mean_K; % average kinetic energy in MeV I.mean_E; % average total energy in MeV I.sigma_x; % mm I.sigma_y; % mm I.sigma_t; % mm/c I.sigma_xp; % mrad I.sigma_yp; % mrad I.sigma_xpx; % mm*mrad I.sigma_ypy; % mm*mrad I.sigma_tpt; % mm/c*permille I.sigma_E; % energy spread in MeV I.sigma_P; % momentum spread in MeV/c I.emitt_x; % mm.mrad, normalised emittance 2.5 Inquiring the bunch properties 33 I.emitt_y; % mm.mrad, normalised emittance I.emitt_z; % mm.permille, normalised emittance (sigma_t*sigma_d) I.emitt_4d; % mm.mrad, 4d normalised emittance I.emitt_6d; % mm.mrad, 6d normalised emittance I.alpha_x; I.alpha_y; I.alpha_z; I.beta_x; % m I.beta_y; % m I.beta_z; % m I.rmax; % mm, largest particle’s xy distance from the origin I.transmission; % the total number of particles in the bunch 2.5.2 Bunch6dT_info If Bis a bunch of type Bunch6dT,B.get_info() returns: I = B.get_info(); I.t; % mm/c I.mean_X; % mm, average H position I.mean_Y; % mm, average V position I.mean_S; % mm, average L position I.mean_Px; % MeV/c, average H momentum I.mean_Py; % MeV/c, average V momentum I.mean_Pz; % MeV/c, average L momentum I.mean_K; % average kinetic energy in MeV I.mean_E; % average total energy in MeV I.sigma_X; % mm I.sigma_Y; % mm I.sigma_Z; % mm I.sigma_Px; % MeV/c I.sigma_Py; % MeV/c I.sigma_Pz; % MeV/c I.sigma_XPx; % mm*MeV/c I.sigma_YPy; % mm*MeV/c I.sigma_ZPz; % mm*MeV/c I.sigma_E; % energy spread in MeV I.emitt_x; % mm.mrad normalised emittance I.emitt_y; % mm.mrad normalised emittance I.emitt_z; % mm.permille, normalised emittance (sigma_z*sigma_pt) I.emitt_4d; % mm.mrad, 4d normalised emittance I.emitt_6d; % mm.mrad, 6d normalised emittance I.alpha_x; I.alpha_y; I.alpha_z; I.beta_x; % m 34 Chapter 2. Beam Models I.beta_y; % m I.beta_z; % m I.rmax; % mm, largest particle’s xy distance from the origin I.transmission; % the total number of particles in the bunch 2.6 Bunch persistency 35 2.6 Bunch persistency Both Bunch6d and Bunch6dT can be loaded and saved on disk and exported in multiple ways. 2.6.1 Saving and loading a beam on disk Two methods can be used to save and load the beam: B.save(filename); B.load(filename); The beam is saved as a binary file. This ensures that the saved and loaded beams are bit-wise identical. Note that the binary format used by each architecture to store double-precision numbers is hardwareand architecture-dependent and may vary from computer to computer. For this reason, a particular architecture may be unable to read a file saved on another hardware architecture because their internal representation of double precision numbers is different. 2.6.2 Exporting as DST file One can save the beam in DST binary format: B.save_as_dst_file(filename, frequency_in_MHz); The RF frequency (expressed in MHz) is required. 2.6.3 Exporting as SDDS file One can save the beam in binary format in SDDS files: B.save_as_sdds_file(filename, description); The string description is optional. 2.6.4 Saving the phase space There are also alternative ways to save a beam as a file. For example, one can extract the phase space using get_phase_space() and then save the phase-space matrix using appropriate Octave or Python commands. 36 Chapter 2. Beam Models 2.7 Spin Polarization Tracking Since version 2.5.0, RF-Track can track spin polarization. Each particle has an associated anomalous magnetic moment G and a spin vector S , a unit 3D vector representing the direction of the particle’s spin polarization in its rest frame. The evolution of the spin vector S follows the Thomas–BMT equation: dS dt =Ω×S, where Ωis the spin precession vector Ω=−q mG+1 γB−Gγ γ+1(β·B)β−G+1 γ+1β×E and •qis the particle charge, •mis the particle mass, •γis the Lorentz factor, •G=g−2 2is the particle anomalous magnetic moment, •β=v/cis the normalized velocity, •B,Eare the magnetic and electric fields in the lab frame (in T and V/m). Notice that RF-Track considers the evolution of the spin vector due to both electric E and magnetic Bfields. 2.7.1 Setting the polarization Once a bunch B0 is created, the spin can be set using one of the two methods: B0.set_polarization (anomalous_magnetic_moment, P ); B0.set_polarization (anomalous_magnetic_moment, P, Sx, Sy, Sz); In these methods, the input arguments are: anomalous_magnetic_moment the anomalous magnetic moment G P the beam polarization, as a number or as a matrix Sx,Sy,Sz the spin polarization vector In the first form, P can be a matrix with three columns and as many rows as particles in the beam. In this case, each row indicates the polarization vector of each particle, and the method sets the spin of each particle to the input. If P is a real scalar number between 0 and 1, then the particles of the bunch are initialized with a degree of vertical polarization equal to P. In the second form, P is the degree of polarization, whereas Sx , Sy , and Sz indicate an arbitrary spin polarization vector. For example, with P=0.8 and Sx=Sy=0 , and Sz=1 , one creates a beam with 80% longitudinal polarization: % Define the parameters A = RF_Track.electron_anomalous_magnetic_moment; % A = (g-2)/2 P = 0.8; % 80% polarization Sx = 0; 2.7 Spin Polarization Tracking 37 Sy = 0; Sz = 1; % Set the polarization of B0 B0.set_polarization (A, P, Sx, Sy, Sz); The anomalous_magnetic_moment is a user-defined quantity. Predefined values are listed in Table 2.3. Variable Name Particle Value RF_Track.electron_anomalous_magnetic_moment electron 0.00115 RF_Track.proton_anomalous_magnetic_moment proton 1.792 RF_Track.muon_anomalous_magnetic_moment muon 0.00116 Table 2.3: Pre-defined anomalous magnetic moments Gof the electron, proton, and muon. Important Remark In RF-Track v2.4.1, the polarization vector of each particle is automatically normalized to 1. Starting with v2.4.2, however, this automatic normalization has been removed. Now, the input polarization of a particle can be any number between 0 and 1. The physical implications of this modification are as follows: • RF-Track v2.4.1 assumes that the spin vector of each macroparticle has a magnitude of 1. This means that each macroparticle is fully polarized in a specific direction. In other words, all of the real particles within the macroparticle have the same identical spin orientation. • RF-Track v2.4.2 and later versions allow the polarization vector of each macroparticle to be less than 1 (still within the range of 0 to 1). This means that the spin vector of each macroparticle is the average polarization of the real particles it represents. Note that any input polarization larger than 1 will still be reduced to 1. 2.7.2 Inquiring the polarization The polarization vector of each particle in a bunch can be accessed using get_phase_space() and the special identifiers %Sx , %Sy , and %Sz . The average polarization of a bunch can be tracked through a beamline using the transport_table with the special identifiers %mean_Sx , %mean_Sy, and %mean_Sz. % Polarization vector of each particle S1 = B1.get_phase_space (’%Sx %Sy %Sz’); % After tracking in a lattice L, % get the average vertical polarization along S T = L.get_transport_table (’%S %mean_Sy’); 2.7.3 Example of spin polarization tracking This example tracks the polarization of a bunch through a solenoid magnet designed to rotate the spin from vertical to horizontal. These solenoid and beam parameters are taken from the design of the Spin Rotator described in the ILC Technical Design Report. 38 Chapter 2. Beam Models 1% Initial phase space, all particles at Z = -5m with momentum 5 GeV/c 2B = zeros(1000, 6); % 1000 rows x 6 columns 3B(:,5) = -5000; % [mm] Z = -5m 4B(:,6) = 5000; % [MeV/c] P = 5 GeV/c 5 6% Create a Bunch6dT for Volume 7B0 = Bunch6dT (RF_Track.electronmass, 0.0, -1, B); 8 9% Set the polarization to 0.8 vertical 10 G = RF_Track.electron_anomalous_magnetic_moment; 11 P = 0.8; 12 Sx = 0; 13 Sy = 1; 14 Sz = 0; 15 B0.set_polarization(G, P, Sx, Sy, Sz); 16 17 % ILC Spin Rotator Solenoid 18 Lsol = 8.32; % m, solenoid length 19 Bsol = 3.15; % T, solenoid on-axis field 20 Rsol = 0.1; % m, solenoid aperture radius 21 Sol = Solenoid (Lsol, Bsol, Rsol); 22 23 % Create a Volume with the ILC Spin Rotator Solenoid 24 V = Volume(); 25 V.set_s0 (-5) % Left boundary 26 V.set_s1 (+5) % Right boundary 27 V.add (Sol, 0, 0, 0, ’center’)% Solenoid centered at 0,0,0 28 29 % Integration parameters 30 V.dt_mm = 10; 31 V.odeint_algorithm = ’rk2’; 32 V.tt_dt_mm = 50; 33 34 % Perform tracking 35 B1 = V.track(B0); 36 37 % Make plots 38 figure(1) 39 hold on 40 S = V.get_transport_table(’%mean_S %mean_Sx %mean_Sy %mean_Sz’); 41 plot(S(:,1) / 1e3, S(:,2), ’displayname’,’S_x’,’linewidth’, 2) 42 plot(S(:,1) / 1e3, S(:,3), ’displayname’,’S_y’,’linewidth’, 2) 43 plot(S(:,1) / 1e3, S(:,4), ’displayname’,’S_z’,’linewidth’, 2) 44 legend(’box’,’off’) 45 xlabel(’S [m]’) 2.7 Spin Polarization Tracking 39 46 ylabel(’Polarization’) The result of this example is visible in Fig. 2.1. Figure 2.1: Example of Spin Tracking. The polarization is moved from vertical to horizontal. 3. Beam Tracking 3.1 Tracking environments RF-Track offers two distinct environments to track particles: the Lattice and the Volume. The environment Lattice represents the accelerator as a list of consecutive elements. Lattice works with Bunch6d and transports the beam, element by element, from the entrance to the exit plane of each element. The environment Volume provides more flexibility than Lattice : it can simulate elements with arbitrary position and orientation in the three-dimensional space and allows elements to overlap. Volume works with Bunch6dT , which is more suitable for space-charge calculations. The possibility for Bunch6dT to handle particle creation at any time and location also allows the simulation of cathodes, field emission, and dark currents. In a Volume , particles can propagate in any direction (even backwards). Several “special” elements, unavailable to Lattice , allow taking full advantage of the flexibility of Volume : e.g., analytic coils and solenoids, where the magnetic field is computed from analytic formulæ and permeates the whole 3D space, allowing for the simulation of realistic fringe fields. It must be mentioned that the environment Lattice can also apply space-charge effects to the beam. However, since Bunch6d maintains the distribution of the particles on the same longitudinal plane, the calculation of the space-charge force needs an on-the-fly extrapolation of each particle’s longitudinal position, using the arrival time and the velocity to reconstruct the three-dimensional spatial distribution. Because of this somehow nonphysical manipulation of the phase space, we deem Lattice more suitable for sections of the accelerator free of space-charge effects. On the contrary, Volume –which through Bunch6dT maintains the full spatial distribution of the beam– should be the preferred choice in space-charge-dominated regimes. These considerations and the capability to superimpose elements make Volume perfect for the simulation of injectors, where a solenoidal magnetic field typically surrounds the gun’s acceleration field, and space charge effects are critical. The additional flexibility of Volume comes at the cost of being computationally more expensive. Summarising, Lattice is straightforward and fast, and its use is recommended in space-charge- 48 Chapter 3. Beam Tracking Table 3.1: A list of all the TrackingOptions odeint_algorithm The integration algorithm. See the following paragraphs for a list of the available algorithms. Default = ‘rk2’ [STRING] odeint_epsabs Absolute error tolerance (’rk’ algorithms only). [REAL] Default = 0.001 odeint_epsrel Relative error tolerance (’rk’ algorithms only). [REAL] Default = 0.0 dt_mm The integration step, the suffix mm stresses the units [mm/c] t_max_mm Sets the max tracking time. Default = +∞[mm/c] t_min_mm Sets the min tracking time (used for backtracking). Default = −∞[mm/c] sc_dt_mm Apply a space-charge kick every sc_dt_mm mm/c[mm/c] cfx_dt_mm Apply the other collective effects every cfx_dt_mm mm/c[mm/c] tt_dt_mm Add an entry of the Transport table every tt_dt_mm mm/c[mm/c] tt_select Select what particles should be included in the transport table computation. It can be one among: ‘ all ’, ‘ active ’, ‘ all_in_volume ’, or ‘ active_in_volume ’. ‘ all ’ also includes particles that haven’t yet come into existence, for instance, particles at the cathode in a photoinjector simulation. [STRING] Default = ‘all’ emission_nsteps If sc_dt_mm is set, applies emission_nsteps space-charge kicks during bunch emission [INTEGER] Default = 10 emission_range Maintain an emission tracking mode for as much time as emission_range times the emission time [REAL] Default = 10.0 wp_dt_mm “wp” stands for “watch point”. Save the beam on disk every wp_dt_mm mm/c [mm/c] wp_basename Base name for the files on disk. Default = “watch_beam”. [STRING] The files will have names: “ watch_beam.XXXXXXX.txt ” and contain 10 columns: 1. Xhorizontal position [mm] 2. Pxhorizontal momentum [MeV/c] 3. Yvertical position [mm] 4. Pyvertical momentum [MeV/c] 5. Zlongitudinal position [mm] 6. Pzlongitudinal momentum [MeV/c] 7. mmass [MeV/c2] 8. Qsingle-particle charge [e] 9. Nnumber of real particles per macroparticle [REAL] 10. particle id [INTEGER] wp_gzip Automatically compress the watch-point files. [BOOLEAN] Default = false verbosity Verbosity level during tracking. 0 = silent; 1 = talkative; 2 = more talkative. Default = 0. [INTEGER] 3.3 Volume 49 Algorithm Description ’analytic’ This algorithm solves analytically the equations of motion assuming a locally constant field. In truly constant fields, this algorithm gives a solution that is exact. In non-constant fields, its accuracy depends on the size of the integration step. This algorithm is symplectic only in truly constant fields. ’leapfrog’ This algorithm updates positions and momenta at different interleaved time points. It’s a second-order symplectic algorithm. ’rk2’ Explicit embedded Runge-Kutta (2, 3) method. ’rk4’ Explicit 4th order (classical) Runge-Kutta. Error estimation is performed using the step-doubling method. ’rkf45’ Explicit embedded Runge-Kutta-Fehlberg (4, 5) method. This method is a good general-purpose integrator. ’rkck’ Explicit embedded Runge-Kutta Cash-Karp (4, 5) method. ’k8pd’ Explicit embedded Runge-Kutta Prince-Dormand (8, 9) method. ’msadams’ A variable-coefficient linear multistep Adams method in Nordsieck form. This stepper uses explicit Adams-Bashforth (predictor) and implicit Adams-Moulton (corrector) methods in P(EC)mfunctional iteration mode. Method order varies dynamically between 1 and 12. The first three are sufficient in most cases; however, we recommend conducting convergence studies and speed tests to select the most suitable algorithm and step size for your specific case. Time limits of the integration In Volume, tracking continues until the absolute time, t , is included between t_min_mm and t_max_mm . The default value for t_max_mm is +∞ , for t_min_mm is −∞ , which means that tracking (or backtracking) continues until all particles have left the Volume (or have reached their creation time, if backtracking). When t_max_mm or t_min_mm is set by the user to a finite number, tracking will stop when treaches it. See the next paragraphs for more details on the stop criterion in Volume. Setting the tracking options A Volume itself is an instance of the structure TrackingOptions ; this allows one to set the volume’s tracking options in this possibly clearer way: % Create Volume V = Volume(); % Set the tracking options V.odeint_algorithm = ’rkf45’; V.dt_mm = 0.1; % mm/c V.verbosity = 2; % Add a few elements V.add (ELEMENT1, 0, 0, 0); V.add (ELEMENT2, 0, 0, 0); V.add (ELEMENT3, 0, 0, 0); 50 Chapter 3. Beam Tracking % Perform tracking V.track(B0); The planes s0and s1 Two planes, called s0 and s1 , orthogonal to the longitudinal axis Z , define the longitudinal extension of a Volume. When one adds elements to a Volume, RF-Track automatically moves s0 backwards and s1forward to accommodate the new element. The planes s0 or s1 can have any position and orientation in space. One can arbitrarily place these planes to achieve specific configurations, or embed a Volume into a Lattice (see next subsection), in multiple ways: % Setting S0 V.set_s0 (z) V.set_s0 (P0, t); V.set_s0 (x, y, z, roll=0, pitch=0, yaw=0); % Setting S1 V.set_s1 (z); V.set_s1 (P0, t); V.set_s1 (x, y, z, roll=0, pitch=0, yaw=0); % Setting one w.r.t. the other V.set_s0_from_s1 (P0, l); V.set_s1_from_s0 (P0, l); In these commands, the arguments are: x, y, z coordinates [m] tthe tracking time [mm/c] la distance [m] P0 a reference particle [Bunch6dT] roll the rotation angle around the Zaxis [rad] pitch the rotation angle around the Xaxis [rad] yaw the rotation angle around the Yaxis [rad] Notice that adding new elements to a Volume automatically updates s1 . So, if you add new elements after set_s1(),s1will likely be changed to include the new elements. End-of-tracking condition Tracking in Volume continues as long as particles exist between these two planes, or until the maximum tracking time (tracking option t_max_mm ) is reached. By default, the maximum tracking time t_max_mm is set to +∞ ( Inf in Octave, or numpy.inf in Python). Beware that if a particle gets trapped in the Volume, tracking will continue indefinitely unless the max tracking time is reached t_max_mm. 3.3 Volume 51 3.3.3 Volume as a Lattice element Once a Volume has been initialized, adding all its elements and setting its planes s0 and s1 , one can insert it into a Lattice and treat it as a standard Lattice element. Figure 3.1 gives a pictorial repreS0 S1 By Lattice 2 Lattice 1 Figure 3.1: A Volume can be sandwiched between two Lattices. sentation of a potential simulation scenario. An illustrative example script follows: % A reference particle placed in (0, 0, 0) with Pz = 100 MeV/c P0 = Bunch6dT (mass, 0.0, +1, [ 0 0 0 0 0 100 ]); % A Volume containing the field map in the figure Dipole = Volume(); Dipole.dt_mm = 1.0; Dipole.odeint_algorithm = ’rk2’; Dipole.add (DIPOLE_MAP, 0, 0, 0, ’center’); Dipole.set_00 (P0, -150); % track backward P0 by 150 mm/c to set s0 Dipole.set_s1 (P0, +150); % track forward P0 by 150 mm/c to set s1 Dipole.set_tt_nsteps(20); % A Lattice L = Lattice(); L.append (Lattice1) L.append (Dipole) L.append (Lattice2) When tracking through L is performed, the particles coming from Lattice1 are distributed over the plane s0 of Dipole , tracked through the Volume Dipole , and collected at s1 to form a Bunch6d suitable to continue tracking through Lattice2. 3.3.4 Transport table and Screens For a detailed tracking of the transport table quantities, a Volume can be sliced into slices using the method set_tt_nsteps() like any Lattice element. Continuing with the example shown in Fig. 3.1, and described conceptually in the previous script, the line 52 Chapter 3. Beam Tracking Dipole.set_tt_nsteps(20); tells RF-Track that we want to track the lattice transport table in 20 steps. To calculate the Lattice transport table, RF-Track will then distribute 20 screens orthogonal to the curved trajectory between s0and s1. 3.4 Particle losses Particles that are lost during tracking can be retrieved from the Lattice or Volume in which they were lost. This can be done using the following methods: M = V.get_lost_particles(); M = L.get_lost_particles(); Both methods return an 11-column matrix with the information on where and when each particle was lost. The information is given in the reference frame of the element itself. In the case of Lattice, the 11 columns are: 1. Xin mm 2. XP in mrad 3. Yin mm 4. YP in mrad 5. Tin mm/c 6. Pin MeV/c 7. Sin m, the longitudinal position at which the particle was lost 8. MASS, in MeV/c2 9. Qin e 10. Nthe macro-particle charge 11. ID the particle ID In the case of Volume, the 11 columns are: 1. Xin mm 2. Px in MeV/c 3. Yin mm 4. Py in MeV/c 5. Zin mm 6. Pz in MeV/c 7. Tin mm/c, the time at which the particle was lost 8. MASS in MeV/c2 9. Qin e 10. Nthe macro-particle charge 11. ID the particle ID 3.5 Synchronization of time-dependent elements with the beam Time-dependent elements -such as RF accelerator structures (whether field maps, travelling-wave, or standing-wave elements), Screen s, and the LaserBeam element used in Compton scattering simulation— require knowledge of the absolute arrival time of the beam in order to operate correctly 3.5 Synchronization of time-dependent elements with the beam 53 and synchronously with the beam. The method set_t0() , available for all time-dependent elements, allows the user to define this arrival time manually. While determining the exact arrival time may seem straightforward for ultra-relativistic beams, where ∆z≈c∆t , it becomes considerably more complex for long beamlines or particles traveling at subluminal velocities ( v<c ), such as heavy ions, protons, or low-energy electrons. To simplify this task, RF-Track provides the autophase() method, common to both Lattice and Volume , which automatically sets the correct synchronization parameters for all time-dependent elements. 3.5.1 Autophasing The autophase() method, available for both Lattice and Volume , synchronizes each timedependent beamline element with the beam by propagating a test bunch through the structure and recording its arrival time, element by element. Below is an example where P0 , a Bunch6d representing an electron bunch with a single reference particle with total charge of 100 pC and average momentum Pref (called Pref in the example), is used for autophasing a lattice L: % Create a reference particle P0 = Bunch6d(electronmass, 100*pC, -1, [ 0 0 0 0 0 Pref ]); % Autophase the lattice Pfinal = L.autophase(P0); In this example, autophase() is given a single reference particle as input. If a full bunch is provided instead, the average particle will be used for autophasing. After calling autophase() , all time-dependent elements are synchronized with the actual arrival time of P0 at each location, and no manual t0 assignment is necessary. Subsequent uses of the lattice L will retain these phases and arrival times. This function returns the final momentum of the reference particle after synchronization. For RF elements, autophase() not only sets the arrival time but also determines the synchronous phase of the structure. Specifically, it adjusts the RF phase so that an input user phase value of phid = 0 (which is the default, if unspecified) corresponds to on-crest acceleration. This feature relieves the user of the tedious task of manually determining the correct RF phase for each structure1. Three important remarks REMARK 1: autophase() should be considered an essential part of the initialization process for aLattice or Volume. It must be called before any tracking is performed. REMARK 2: Any time-dependent element whose reference time is set using set_t0() before being added to a Lattice or Volume will be excluded from the autophasing process and will retain its manually defined reference time. REMARK 3: If the user does not explicitly call autophase() , RF-Track will automatically invoke it during the first call to the track() method. 1In fact, standard electromagnetic simulation codes provide field maps with an arbitrary RF phase. 54 Chapter 3. Beam Tracking 3.5.2 Automatic setting of magnetic strengths In RF-Track, magnetic fields are typically defined using either the magnetic strength or the gradient. This is the case, for instance, with the Quadrupole and Multipole elements discussed in the next chapter. However, when designing accelerator structures such as linacs composed of FODO cells interleaved with RF cavities, it is often challenging to specify the exact magnet gradients. This is because the magnetic strength depends on the beam momentum, which typically varies along the linac and is not known a priori. In such situations, it is more practical to define the magnets using their normalised strengths — quantities that are independent of energy and determined solely by the optical properties of the lattice. For example, the quadrupole’s parameter k1 [1/m −2 ] describes the normalised strength of a quadrupole. When k1 is provided, RF-Track can defer the computation of the actual magnetic gradient until the beam energy is known. This is achieved by setting the reference momentum to NaN when creating the magnet. For example: % We know the normalized focusing strength, k1, % but the beam rigidity is not yet defined. k1 = 0.12; % m^-2, normalised strength L = 0.1; % m, quadrupole length P_Q = NaN;% Beam momentum is unknown at this stage % Define a quadrupole with deferred gradient computation Q = Quadrupole (L, P_Q, k1); The quadrupole gradient will then be automatically computed and assigned during the call to autophase() . If the quadrupole length is set to zero, then k1 is interpreted as the integrated strength k1L[1/m], allowing the simulation of thin-lens quadrupoles. 3.6 Backtracking RF-Track supports beam backtracking. During backtracking, element misalignments and all deterministic collective effects are fully accounted for. This includes effects such as space charge or wakefields, but excludes stochastic processes like multiple Coulomb scattering or quantum incoherent synchrotron radiation emission. Backtracking can be performed using the btrack() method, which is available for both Lattice and Volume objects. For example: % B1 is the desired final phase-space distribution B1 = Bunch6d(...); % Backtrack to find B0, the initial distribution that evolves into B1 B0 = L.btrack(B1); This feature is particularly useful when the goal is to determine the initial beam distribution that yields a specified final distribution at the end of a beamline. 4. Beamline Elements 4.1 Introduction RF-Track offers a comprehensive set of elements, ranging from conventional matrix-based quadrupoles and sector bends to more sophisticated options, such as complex field maps and analytic fields that permeate the entire 3D space. This chapter gives a short description of the constructors for each element. In some elements, specific attributes cannot be directly set via the constructor but can be set using appropriate “set” methods. The most relevant methods are given for each element type. The user can access a list of all methods using the interactive prompts of Octave and Python. 4.1.1 Methods available to all elements Many practical methods are available to all elements to achieve specific setups. Setting the element name The element can have an optional name that can be used to look up the elements in the tracking environments. Two methods allow to set and get an element’s name: E.set_name (STRING); E.get_name (); Setting the aperture All elements are equipped with an aperture, which, by default, is not considered. The aperture can be set and inquired using the following set of methods: E.set_aperture_x (Rx); % m E.set_aperture_y (Ry); % m E.set_aperture_shape (SHAPE); % m E.set_aperture (Rx, Ry, SHAPE); % m E.get_aperture_x (); % return m 56 Chapter 4. Beamline Elements E.get_aperture_y (); % return m E.get_aperture_shape (); % return a string The arguments are: Rx,Ry the aperture radii [m] SHAPE a string among ’none’, ’rectangular’, and ‘circular’ [STRING] The ‘circular’ aperture becomes elliptical when Rx 6=Ry . The aperture is checked at each integration step. Adding constant magnetic and electric fields Available to most element types, two methods allow embedding an element in a constant magnetic or electric field with any orientation in space: E.set_static_Bfield (Bx, By, Bz); % T E.set_static_Efield (Ex, Ey, Ez); % V/m Where Eis the instance of an element. Activating collective effects Collective effects can be added to any element using the following method: E.add_collective_effect (CFX); See the dedicated chapter for a list of the collective effects implemented in RF-Track. 4.1.2 Inquiring and plotting the electromagnetic field Each RF-Track element can be queried for the electromagnetic field it represents at any point. The method get_field() can be used to do this. See the following example: Q = Quadrupole (1.0, 100.0); % Or any other element % Inquire the % field at point (x,y,z,t) x = 0; % mm y = 10; % mm z = 20; % mm t = 30; % mm/c [E,B] = Q.get_field (x, y, z, t); % Returns E in [V/m] and B in [T]. The input arguments of get_field() can be either scalars or one-dimensional vectors to query a sequence of points at once. E and B are three-element vectors containing Ex , Ey , Ez and Bx , By , Bz respectively in the case of a single point, or three-column matrices with as many rows as there are points queried. Compound elements such as Lattice and Volume can also be queried for the field. This makes it possible, for example, to plot the total electromagnetic field of all elements in a given environment at any point. 4.1.3 Tracking the beam properties Both Lattice() and Volume() offer the possibility of storing average beam quantities, such as beam size, emittance, energy spread, dispersion, etc., during tracking into a “transport table”. After 4.1 Introduction 57 tracking, such a transport table can be retrieved using the method get_transport_table() available to both Lattice and Volume. Like get_phase_space() for the beam, get_transport_table() allows the user to inquire about specific quantities. In Lattice() , the user can choose the number of points at which the phase space is sampled using the element’s method set_tt_nsteps(). For example: D = Drift (1); % 1-meter-long drift D.set_tt_nsteps (100); which will track the beam through the drift D in 100 steps, sampling the phase space in that many points. Table 4.1 lists all the accepted identifiers. To enable a transport table in Volume() , it is sufficient to specify the option tt_dt_mm in the tracking options, specifying the time interval, in mm/ c , between two consecutive samplings. Table 4.2 lists all the accepted identifiers. 64 Chapter 4. Beamline Elements Set methods The obvious ones, plus: S.set_h(H); % 1/m, set the curvature of the reference system S.set_hgap(HGAP); % m, the half gap of the magnet S.set_fint(FINT); % the fringe field integral The parameter FINT follows the same definitions as in MAD-X: FINT =Z∞ −∞ By(s)(B0−By(s)) gB2 0 ds, with g=2HGAP. The default value FINT of zero corresponds to the hard-edge approximation, i.e. a rectangular field distribution. For other approximations, one can refer to the following pre-computed values of FINT: Linear field drop-off 1/6 Clamped “Rogowski” fringing field 0.4 Unclamped “Rogowski” fringing field 0.7 “Square-edged” non-saturating magnet 0.45 4.2 Matrix-based elements 65 4.2.4 Rectangular bending dipole RBend is a rectangular bend, a bending magnet whose entrance and exit windows are parallel. Its reference system is curved by default, and its length is intended to be the length of the straight line joining the entry and exit points. Constructor R = RBend (L, angle, P_Q, E1=0, E2=0); The arguments are: Lthe rectangular bending magnet length [m] angle the bending angle [rad] P_Qthe beam’s magnetic rigidity, P/q[MV/c] E1,E2 the entrance and the exit angles [rad] As the RBend is internally implemented as a SBend with appropriate entrance and exit edge angles, this element features the same “Set” and “Get” methods of a sector bend. 66 Chapter 4. Beamline Elements 4.2.5 Corrector The element Corrector creates a magnetic steerer for orbit correction. Constructors C = Corrector (L); C = Corrector (L, Kx, Ky); The input arguments describing the Corrector are: Lthe corrector length [m] Kx,Ky the integrated horizontal and vertical corrector strengths [T.mm] 4.3 Special elements 67 4.3 Special elements In RF-Track, special elements represent components or devices that cannot be implemented using simple linear transfer matrices, but require numerical integration in non-linear fields. These include coils, solenoids, various RF fields, undulators, etc. 4.3.1 Coil The element Coil creates the magnetic field generated by an electromagnetic coil. Constructors C = Coil (L, B0, R ); The input arguments describing the coil are: Lthe element length [m] B0 the peak on-axis field [T] Rthe coil radius. [m] The coil is placed in the middle of the specified length L . When a Coil is placed in a Volume, its field permeates the whole 3D space. Warning: When used in Lattice, the field exists only within the extent of the specified element length, with the coil placed in the middle. See also the element Solenoid. 68 Chapter 4. Beamline Elements 4.3.2 Solenoid The element Solenoid allows the insertion of a Solenoid magnet into a Lattice or a Volume. When a Solenoid is inserted into a Lattice, RF-Track will treat it as a matrix-based element. When Solenoid is inserted into a Volume, RF-Track computes the 3D magnetic field of the magnet using analytic equations. In this case, the field includes realistic fringe regions and permeates the whole 3D space, overlapping with other fields in the same Volume. Constructor S = Solenoid (L=0, B0=0, R=0); S = Solenoid (L, B0, Rmin, Rmax, nsheets); The arguments are: Lthe length of the solenoid [m] B0 the on-axis peak field, [T] Rthe solenoid radius [m] Rmin inner radius [m] Rmax outer radius [m] nsheets number of current sheets [INTEGER] In both Lattice and Volume, the argument radius, R , sets the aperture radius of the magnet. In Volume, its value is also used to compute the 3D field. Figure 4.2 shows the Solenoid field computed by RF-Track. 1.00 0.75 0.50 0.25 0.00 0.25 0.50 0.75 1.00 z [m] 0.6 0.4 0.2 0.0 0.2 0.4 0.6 x [m] 0.5 1.0 1.5 2.0 2.5 | B | [T] Figure 4.2: The Solenoid field computed by RF-Track. 4.3 Special elements 69 4.3.3 Undulator The element Undulator allows the insertion of a planar Undulator into a Lattice or a Volume. When particles travel through an Undulator , RF-Track computes the 3D magnetic field of the magnet using analytic equations and integrates the equations of motion using numerical integration. Constructor U = Undulator (lperiod, K, nperiods, kx2=0); The arguments are: lperiod the length of an undulator period [m] Kthe undulator K parameter in xdirection [-] nperiods the number of periods [INTEGER] kx2 the curvature of the pole surface, k2 x(default 0) [1/m2] When kx2 >0 , the poles are bent inwards; when kx2 <0 , the poles are bent outwards. See Figure 4.3. <latexit sha1_base64="gj0BCCXOgkMQsH16hQ6f9e47E3k=">AAAB8XicbVBNS8NAEJ3Ur1q/oh69LBbBU0mKVE9S9OKxgv3ANpbNdtMu3WzC7kYsof/CiwdFvPpvvPlv3LQ5aOuDgcd7M8zM82POlHacb6uwsrq2vlHcLG1t7+zu2fsHLRUlktAmiXgkOz5WlDNBm5ppTjuxpDj0OW374+vMbz9SqVgk7vQkpl6Ih4IFjGBtpPtx/+mhii6RU+rbZafizICWiZuTMuRo9O2v3iAiSUiFJhwr1XWdWHsplpoRTqelXqJojMkYD2nXUIFDqrx0dvEUnRhlgIJImhIazdTfEykOlZqEvukMsR6pRS8T//O6iQ4uvJSJONFUkPmiIOFIRyh7Hw2YpETziSGYSGZuRWSEJSbahJSF4C6+vExa1Ypbq9Ruz8r1qzyOIhzBMZyCC+dQhxtoQBMICHiGV3izlPVivVsf89aClc8cwh9Ynz9ohY9z</latexit> k2 x>0 <latexit sha1_base64="5oe3VtZEkdrJSQAyK7+S5M+EkpA=">AAAB8XicbVBNS8NAEJ3Ur1q/oh69LBbBU0mKVA8eil48VrAf2May2W7apZtN2N2IJfRfePGgiFf/jTf/jZs2B219MPB4b4aZeX7MmdKO820VVlbX1jeKm6Wt7Z3dPXv/oKWiRBLaJBGPZMfHinImaFMzzWknlhSHPqdtf3yd+e1HKhWLxJ2exNQL8VCwgBGsjXQ/7j89VNElckp9u+xUnBnQMnFzUoYcjb791RtEJAmp0IRjpbquE2svxVIzwum01EsUjTEZ4yHtGipwSJWXzi6eohOjDFAQSVNCo5n6eyLFoVKT0DedIdYjtehl4n9eN9HBhZcyESeaCjJfFCQc6Qhl76MBk5RoPjEEE8nMrYiMsMREm5CyENzFl5dJq1pxa5Xa7Vm5fpXHUYQjOIZTcOEc6nADDWgCAQHP8ApvlrJerHfrY95asPKZQ/gD6/MHZXePcQ==</latexit> k2 x<0 Figure 4.3: The parameter kx2 describes the shape of the undulator’s poles. 70 Chapter 4. Beamline Elements 4.3.4 Transfer Line The TransferLine element can be used to transport a beam through an entire beamline without having to specify each element individually, only by giving pairs of Twiss parameters. In the transverse plane, the element tracks the beam between each pair of Twiss parameters {β,α}1 and {β,α}2using a transfer matrix in Twiss form: M1→2= qβ2 β1(cos µ+α1sin µ)pβ1β2sin µ (α1−α2)cos µ−(1+α1α2)sin µ √β1β2qβ1 β2(cos µ−α2sin µ)  with phase advance µ . If horizontal and vertical chromaticities are provided, the phase advance in both planes is adjusted accordingly. In the longitudinal plane, a drift is applied with length L=L0(1+αCδ) , where L0 is the length of the element’s section, αC is the momentum compaction factor, and δis the relative momentum difference. This element has three modes of operation: • If one gives two sets of Twiss parameters, RF-Track will use a Twiss matrix. • If one gives a Twiss table or a Twiss file from MAD-X, RF-Track will track through each consecutive pair of Twiss parameters using the Twiss transfer matrix. • If one gives just one set of Twiss parameters, RF-Track will match any bunch being tracked to the desired set of Twiss parameters. Like any other Lattice elements, a TransferLine can be segmented into steps, and collective effects can be applied through the TransferLine. Constructors T = TransferLine ("twiss_file.dat", Pref); T = TransferLine (twiss_matrix, Pref, DQx, DQy, momentum_compaction); T = TransferLine (twiss_matrix, Pref); The arguments are: "twiss_file.dat" a Twiss file from MAD-X twiss_matrix a matrix containing the Twiss parameters. This matrix can be either a 7-column matrix, where each row contains: [S,βx,αx,µx,βy,αy,µy] or an 11-column matrix, where each row contains: [S,βx,αx,µx,βy,αy,µy,Dx,Dpx,Dy,Dpy ] S,βx,βy,Dx,Dyare expressed in [m] µx,µyare expressed in [2π] Pref the reference momentum [MeV/c] DQx,DQy the horizontal and vertical chromaticities [2π] momentum_compaction the momentum compaction All the input quantities follow their respective MAD-X definitions. Note: This element works in Lattice only. 4.3 Special elements 71 Example of matching section We propose the example of a matching section implemented as a TransferLine , matching a set of Twiss parameters [β0,x,α0,x,β0,y,α0,y] into [β1,x,α1,x,β1,y,α1,y] , with phase advance µx and µyover a length L: % Input data Pref = 100; % MeV/c, reference momentum L=5;% m, length of the matching section % Initial Twiss parameters b0_x = 1; % m, horizontal beta function b0_y = 2; % m, vertical beta function a0_x = a0_y = 0; % alpha % Final Twiss parameters b1_x = 5; % m, horizontal beta function b1_y = 10; % m, vertical beta function a1_x = a1_y = 0; % alpha % Phase advance between initial and final mu_x = pi/2; % phase advance in x mu_y = pi/3; % phase advance in y % Input 2-row matrix for the constructor T=[0,b0_x, a0_x, 0.0, b0_x, a0_y, 0.0; % 1st row, initial Twiss L, b1_x, a1_x, mu_x, b1_x, a1_y, mu_y];% 2nd row, final Twiss % Transfer matrix for Lattice M = TransferLine(T, Pref); 72 Chapter 4. Beamline Elements 4.3.5 Travelling-wave structure The element TW_Structure allows the simulation of the TM01n modes in a travelling-wave (TW) structure using an analytic description of the field, as in the equations below. A TW structure is a metallic enclosure where electromagnetic fields travelling with a specific phase advance along the structure can be excited. Figure 4.4 shows a model example of a travelling-wave structure. Figure 4.4: A model of a metallic travelling-wave structure. Table 4.3: Fourier series expansions of the TM01nelectromagnetic fields in a TW structure Quantity Value k0 ∆ϕ L knk0+2nπ L βn ω ckn qnsω c2−(kn)2 Ez(r,z,t) ∞ ∑ n=−∞ ansin[(ωt+φ0)−knz]×(J0(qnr)|βn|≥1 I0(qnr)|βn|<1 Er(r,z,t) ∞ ∑ n=−∞ ankn qn cos[(ωt+φ0)−knz]×(J1(qnr)|βn|≥1 I1(qnr)|βn|<1 Bθ(r,z,t) ∞ ∑ n=−∞ anq2 n+k2 n ωqn cos[(ωt+φ0)−knz]×(J1(qnr)|βn|≥1 I1(qnr)|βn|<1 For such structures, if L and ∆ϕ indicate the length and phase advance of one cavity cell, then the Fourier series expansions of the excited TM01n modes in the cell are given by the expressions presented in Table 4.3. Constructors T = TW_Structure ([ an, ... ], n, freq, ph_advance, number_of_cells); T = TW_Structure (); The arguments are: 4.3 Special elements 73 [an, ... ] a vector with the Fourier coefficients (see Table 4.3) [V/m] nthe index of the first coefficient [INTEGER] freq the rf frequency [Hz] ph_advance the structure’s phase advance per cell [rad] number_of_cells the number of cells in the structure. [REAL] A positive number of cells indicates that the structure starts from the middle of the cell. A negative number indicates that the structure starts from the beginning of a cell. Follows an example with three cells: •number_of_cells = 3 •number_of_cells = -3 Example of travelling-wave structure We present an example of an ideal travelling-wave structure with three cells having a phase advance of 2π/3, operating at 3 GHz with a gradient of 10 MV/m. a0 = 10e6; % V/m, gradient freq = 3e9; % Hz ph_adv = 2*pi/3; % rad, phase advance n_cells = 3; % number of cells n=0;% TM 01n mode TW = TW_Structure(a0, n, freq, ph_adv, n_cells); Figure 4.5 presents the longitudinal field profile when one changes the parameter n from 0 to 3 (top to bottom). 80 Chapter 4. Beamline Elements Set methods The obvious ones, plus: S.set_Bn([ B0 B1 B2 ... ]); % T/m^n, set the field derivatives S.set_KnL(P_over_Q, [ K0L K1L K2L ... ]); % set the integrated coefficients S.set_strengths([ S0 S1 S2 ... ]); % MV/c/m^n, set the integrated strengths Two methods allow setting the tracking options: S.set_nsteps(N); % set the number of integration steps [default 10] S.set_odeint_algorithm(name); % e.g. ’rk2’, ’leapfrog’, etc. 4.3 Special elements 81 4.3.9 Absorber The Absorber element is a block of matter where three collective effects affecting the beam are applied simultaneously: Multiple Coulomb scattering, Energy Straggling, and Stopping Power 1 . These enable the simulation of interactions between bunch particles and materials. Optionally, Absorber s can have a circular or a rectangular shape. This is specified using the method: A.set_shape (shape, ax [m], ay [m] ) where shape can be either “ circular ”, or “ rectangular ”; “ ax ” and “ ay ” are the radii in x and y. Note that, when an Absorber is placed in a Volume, the computation of collective effects must be activated for it to be effective. This can be done by setting the TrackingOption cfx_dt_mm , see chapter 3. When an Absorber is placed in a Lattice, these collective effects are enabled by default. Constructors A = Absorber (L, material_name); A = Absorber (L, X0, Z, A, density, I=-1); The input arguments describing the material are: Labsorber length [m] material_name one among: ‘ air ’, ‘ water ’, ‘ beryllium ’, ‘lithium’, ‘liquid_hydrogen’ [STRING] X0 radiation length [cm] Zatomic number of absorber [INTEGER] Aatomic mass of absorber [g mol−1] density density of the absorber [g cm−3] I mean excitation energy. If unspecified or -1, an empirical approximation is used. [eV] Main methods A set of ‘enable’ and ‘disable’ methods allows the full customization of this element. A.enable_log_term(); A.enable_fruehwirth_model(); % DEFAULT A.enable_wentzel_model(); % DEFAULT A.enable_stopping_power(); % DEFAULT A.enable_energy_straggling(); % DEFAULT A.enable_multiple_coulomb_scattering(); % DEFAULT A.disable_log_term(); % DEFAULT A.disable_fruehwirth_model(); A.disable_wentzel_model(); A.disable_stopping_power(); A.disable_energy_straggling(); A.disable_multiple_coulomb_scattering(); 1For a detailed description of the three effects, see the dedicated chapter. 82 Chapter 4. Beamline Elements 4.3.10 Electron cooler RF-Track can simulate electron cooling using a sophisticated hybrid kinetic model in which the beam’s macro-particles interact with a cold, magnetised electron plasma. The electron plasma is modelled with 3D meshes, allowing the user to specify arbitrary transverse and longitudinal variations of density and current. The following panel presents its constructor and the main methods to customize the electron cooler simulation fully: EC = ElectronCooler (L, rx, ry, density, Vz); % Plasma mesh EC.set_temperature (Tr, Tl); EC.set_electron_mesh (Nx, Ny, Nz, density, Vx, Vy, Vz); % uniform plasma EC.set_electron_mesh (Nz, DENSITY2D, VX2d, VY2D, VZ2D); % 2D profile EC.set_electron_mesh (DENSITY3D, VX3D, VY3D, VZ3D); % full 3D mesh % Magnetic field EC.set_static_Bfield (Bx, By, Bz); % Plasma particle (default: electron) EC.set_Q (Q=-1); EC.set_mass (mass=electronmass); The main input arguments are: Lthe electron cooler’s length [m] rx,ry the horizontal and vertical radius of the plasma [m] Tr,Tl the radial and longitudinal temperature of the plasma [eV] density the plasma density [m−3] DENSITY2D the plasma density as a Nx×Ny matrix [m−3] DENSITY3D the plasma density as a Nx×Ny×Nz matrix [m−3] Bx,By,Bz the constant magnetic field magnetizing the plasma [T] Vx,Vy,Vz the velocity of the plasma beam [c] VX2D,VY2D,VZ2D the velocity of the plasma beam as a Nx×Ny matrix [c] VX3D,VY3D,VZ3D the velocity of the plasma beam as a Nx×Ny×Nz matrix [c] 4.3 Special elements 83 4.3.11 Adiabatic matching device An AdiabaticMatchingDevice (AMD) is a magnetic element providing a strong tapered solenoid field, typically used in positron sources. This device, also known as a “Flux Concentrator”, is essential to capture the positrons right after their creation. In RF-Track, this element is implemented as an analytic 3D field. The on-axis field provided by an AMD is purely longitudinal and is described as, Bz(z) = B0 1+αz. Figure 4.7 shows the on-axis field Bz as a function of z for an AMD with B0=7T and α=60m−1 . Constructor AMD = AdiabaticMatchingDevice (L, B0, ALPHA); 0 1 2 3 4 5 6 7 0 0.05 0.1 0.15 0.2 0.25 Bz (T) Z (m) Figure 4.7: Magnetic profile in an AMD. The arguments are: Lthe AMD length [m] B0 the peak on-axis field, B0[T] ALPHA the parameter α[m−1] Main methods The aperture of an Adiabatic Matching Device has the shape of a truncated cone. You can specify the aperture radius of each end of the AMD using the following two methods: A.set_entrance_aperture (R1); A.set_exit_aperture (R2); Where R1 and R2 are the entrance and exit aperture radii in metres. 84 Chapter 4. Beamline Elements 4.3.12 Space-charge Field SpaceCharge_Field is a special element that allows you to know the electromagnetic field generated by an arbitrary distribution of particles at any point in space, using the method get_field() . It can be used to simulate weak-strong interactions in Volume. In its constructor, SpaceCharge_Field accepts a particle distribution given as Bunch6dT as an input. Inside the rectangular bounding box enclosing all particles, the field is computed using the same 3D space-charge calculation routines used during tracking: a 3D PIC based on FFT-integrated retarded Green’s functions in free space. Outside the bounding box (that is, at large distances), the field is computed using a Cartesian multipole expansion of the charge and current distribution in space to 5th order. Since SpaceCharge_Field creates an electromagnetic field that permeates the entire space, this element is intended for use in Volume only. Figure 4.8 provides a visual example of SpaceCharge_Field and the distribution used for its initialisation. Constructor SC = SpaceCharge_Field (B0T, Nx, Ny, Nz, Vz_slices=1); The arguments are: B0T an arbitrary particle distribution [Bunch6dT] Nx, Ny, Nz the number of mesh points for the 3D PIC solver [INTEGER] Vz_slices the number of velocity slices accounting for relativistic effects [INTEGER] Figure 4.8: Electromagnetic field generated by a particle distribution, computed by the element SpaceCharge_Field . The field extends both internally and externally to the particle distribution generating it. A SpaceCharge_Field can be superimposed on external fields for the accurate simulation of sources. (Figure courtesy of Manon Boucard, PMB-Alcen, France). 4.3 Special elements 85 4.3.13 Travelling-Wave Field The element TW_Field allows the insertion of a travelling-wave structure directly from its shunt impedance, group velocity, and quality factor. With this element, explicit field-mapping fields for the electric and magnetic fields aren’t necessary. In order to compute the field, RF-Track also needs to know the input power into the structure. The input power can be provided in two ways: 1. Constant Pin in W. When given a Pin, RF-Track computes the steady state of the field. 2. Dynamic Pin . The input power is provided as a 1D array of time-dependent values, and an injection time ttextin j is specified, determining the moment the bunch enters the structure along the power input. Constructors TW = TW_Field (P_in, Q, r_Q, VG, freq, ph_adv, n_cells, z0_L=0 ); TW = TW_Field (P_in, dt, t_inj, Q, r_Q, VG, freq, ph_adv, n_cells ); The arguments are: P_in The input power in the structure. In the first form of the constructor, it is a scalar number. The assumption here is an infinitely long input power pulse. The field is computed after the filling is completed. In the second constructor, it is a vector of numbers sampling the input power pulse at regular time interval dt [W] dt The time step of the input vector P_in [mm/c] t_inj the time of injection of the beam in the structure, with respect to the input power vector P_in [mm/c] QThe quality factor along the structure r_QNormalised shunt impedance per unit length array (r/Q) [Ω/m] vg Group-velocity array [c] ph_adv The phase advance per cell [rad] n_cells The number of cells [INTEGER] z0_L The length of the input and output couplers, expressed as a fraction of the cell length Lcell . If z0_L6=0 , RF-Track will add a standing-wave coupler before and one after the travelling-wave body of the structure, with a length Lcoupler =z0_L·Lcell, where Lcell is the length of the structure’s cell [REAL] Notice that the key vectors describing the structure field, Q , r_Q , VG , must have the same number of elements. 86 Chapter 4. Beamline Elements 4.3.14 LaserBeam In Lattice, RF-Track can simulate inverse Compton scattering (ICS) between any charged particle and a laser beam in an element called LaserBeam . Chapter 5 illustrates this element and the ICS simulation in detail. 4.3 Special elements 87 4.3.15 Volume as a Lattice element The tracking environment Volume can be inserted into a Lattice like any standard element. When a Bunch6d , which is being tracked through a Lattice, enters a Volume that has been inserted into a Lattice as an element, its particles are placed over the are emitted one by one from the Volume’s surface S0 , according to their time coordinate. Tracking is then performed in the Volume by time integration and stops when no more particles exist between S0and S1. During tracking, when a particle reaches S1 and leaves the Volume, its position at S1 and its arrival time are recorded, and a bunch of type Bunch6d is formed for propagation through the rest of the lattice. This way, a section where time integration is used with all the flexibility offered by Volume can be placed in a Lattice. Note that the planes S0 and S1 can have any orientation in space. This allows for changes in the reference system, for example, when tracking in the field map of a bending magnet or when joining regions at an angle. 88 Chapter 4. Beamline Elements 4.4 Field maps One of RF-Track’s strengths is tracking in field maps. Real and complex field maps, in one, two, or three dimensions, are accepted to simulate static fields, backward— or forward-travelling fields, or standing-wave radiofrequency fields. Elements based on a one-dimensional on-axis field map still provide three-dimensional fields. By proving the longitudinal component of the field along the axis of symmetry, RF-Track extends the field off-axis according to Maxwell’s equations, assuming cylindrical symmetry. Two-dimensional field maps allow the simulation of cylindrical-symmetric fields from the field over a plane. Threedimensional field maps allow tracking in most generic electromagnetic fields. All field maps accept 3D Cartesian mesh grids with regular spacing. Interpolation methods RF-Track offers two interpolation methods • Linear interpolation (LINT), the field is evaluated at any point by linearly interpolating with the eight nearest mesh points. This method is very fast and is the default method. • Cubing interpolation (CINT), the field is evaluated in any point by cubically interpolating with the nearest 64 mesh points. This method is significantly slower but produces a smoother field. Particle losses in field maps In RF-Track, three-dimensional field maps can contain more information than the field itself. They can contain information about the walls of the element in which the field is embedded, using the special value “not-a-number” (NaN). In floating-point arithmetic, NaN is a numeric data type that can be interpreted as an undefined value. RF-Track interprets the presence of NaNs in a field map as a “wall”. If a particle hits a NaN, RF-Track flags the particle as lost. This allows precise detection of losses in 3D space, even in complex geometries. Interpolation method and losses detection In the case of LINT, the interpolation uses the eight closest vertexes of the 3D mesh cell enclosing the point of interest. Notice that, in this case, the granularity of the loss detection coincides with the mesh cell size. In the case of CINT, the interpolation uses the 64 closest vertexes, as it considers the cube with 3×3×3 mesh cells surrounding the point of interest. This means that the granularity of the loss detection is effectively three times the size of a mesh cell. CINT is generally smoother than LINT, as it considers four adjacent points instead of 2. The price for the increased smoothness is computational time. 4.4 Field maps 89 4.4.1 1D RF field maps An RF_FieldMap_1d provides an RF oscillating electromagnetic field based on the on-axis electric field. RF-Track assumes cylindrical symmetry around the structure’s axis and reconstructs the off-axis transverse and longitudinal electromagnetic field, fulfilling Maxwell’s equations. The user must provide a regular Cartesian 1D mesh of the longitudinal electric field Ez as complex (travelling waves) or real (standing waves) numbers. Constructors RF = RF_FieldMap_1d (Ez, hz, length, frequency, direction, P_max = 1, P_actual = 1); RF = RF_FieldMap_1d_CINT ( " ); The default element RF_FieldMap_1d uses linear interpolation along the longitudinal axis and linear extrapolation in the radial direction. The variant, RF_FieldMap_1d_CINT , uses cubic interpolation along the longitudinal axis and cubic extrapolation in the radial direction. The required input arguments are: Ez The on-axis electric field [V/m] hz The mesh cell in the longitudinal direction [m] length the total length of the element (if -1 take the field map length) [m] frequency The RF frequency (use 0 for static fields) [Hz] direction 0 : static field 1 : forward-travelling field −1 : backward-travelling field NOTE: standing-wave fields can have +1 or -1 identically P_map The input power used to generate the input field map (default 1) [W] P_actual The actual input power to operate the element (default 1) [W] If the user only needs to simulate an electric or magnetic field, the constant number zero, 0 , can be provided for the unnecessary field components. Main methods A set of dedicated methods allows the user to set the phase, the reference time, and the actual input power to the structure: RF.set_t0(T0); RF.unset_t0(); RF.set_phid(PHID); RF.set_P_actual(P); RF.set_smooth(N); The input parameters are: 96 Chapter 4. Beamline Elements 4.4.5 2D static magnetic field maps An element of type Static_Magnetic_FieldMap_2d provides a magnetic field based on a field map defined over a 2D plane. RF-Track assumes cylindrical symmetry around the element’s axis. The user must provide a regular Cartesian 2D mesh per each field component: Br and Bz . The two input 2D meshes, Bri j ,Bzi j have two indexes, i j : the first, i , runs along the longitudinal direction, and the second, j , along the radial direction. Constructors S = Static_Magnetic_FieldMap_2d (Br, Bz, hr, hr, length=-1); S = Static_Magnetic_FieldMap_2d_CINT ( " ); The element Static_Magnetic_FieldMap_2d uses linear interpolation. The variant, Static_Magnetic_FieldMap_2d_CINT, provides cubic interpolation. The required input arguments are: Br,Bz 2D mesh of the radial and longitudinal magnetic field [T] hr,hz The sizes of a mesh cell in the radial and longitudinal directions [m] length the total length of the element (if -1 take the field map length) [m] 4.4 Field maps 97 4.4.6 3D static magnetic field maps An element of type Static_Magnetic_FieldMap provides a magnetic field based on a 3D field map. It accepts three input 3D meshes or four if one specifies the field using its scalar and vector potentials. Constructors S = Static_Magnetic_FieldMap (Bx, By, Bz, x0, y0, hx, hy, hz, length); S = Static_Magnetic_FieldMap (Ax, Ay, Az, PhiM, x0, y0, hx, hy, hz, length); The required input arguments are: Bx,By,Bz 3D mesh of the magnetic field [T] Ax,Ay,Az 3D mesh of the vector magnetic potential 0[T m] PhiM 3D mesh of the scalar magnetic potential 0[T m] x0,y0 The position of the bottom-left starting point of the field map in the x−yplane [m] hx,hy,hz The sizes of a mesh cell in the x,y, and zdirections [m] length the total length of the element (if -1 take the field map length) [m] 98 Chapter 4. Beamline Elements 4.5 Beam diagnostics 4.5.1 Beam position monitor The Bpm element adds a beam position monitor to a Lattice. The Bpm can be thick and have a user-defined resolution. The scaling error can also be applied. Reading a Bpm can be done after tracking a beam. Note that subsequent readings for the same Bpms will result in different readings, as each measurement is affected by the resolution error. The Bpm element is currently only implemented in Lattice. Constructor B = Bpm (L, resolution ); The input arguments describing the material are: L the element length. The reading occurs in the middle of the specified length. [m] resolution the resolution [mm] After tracking, to get all Bpms’ readings at once, one can call the method get_bpm_readings() of the Lattice object. Get methods [X,Y] = B.get_reading(); % return mm res = B.get_resolution(); % return mm Set methods B.set_resolution(res); % accept mm B.set_scaling_factor(X_scaling, Y_scaling = X_scaling); 4.5 Beam diagnostics 99 4.5.2 Screens The Screen element adds a Screen to a Lattice or Volume. A screen is a thin element that captures a snapshot of the phase space of a Beam or just a Bunch6d when this is traversing the screen itself. Screens can have a finite extension and a specific time window; see the methods below. When a time window is specified, the screen is automatically synchronized to the first bunch traversing the screen unless the user selects a specific activation time. In a Lattice, Screens can be placed between elements with any arbitrary offset. In Volume, Screens can be placed in any position with any orientation in space. When a bunch traverses a screen, RF-Track retains the arrival time of each particle at the screen and creates a Bunch6d with the bunch’s phase space in the screen’s reference system. After tracking, the phase space of the bunch (or of the beam) can be retrieved from Volume and Lattice using the methods get_bunch_at_screens() or get_beam_at_screens() , depending on whether one is tracking a single or multi-bunch beam. These methods return a list of Bunch6d’s or Beams objects, one per screen. Constructor S = Screen (); Get methods B = S.get_bunch(); % return a Bunch6d B = S.get_beam(); % return a Beam t0 = S.get_t0(); % returns the reference time in mm/c Set methods S.set_width (W); % mm, screen width W S.set_height (H); % mm, screen height H S.set_time_window (T); % mm/c, screen time window T S.set_t0 (t0); % set reference time [default: unset] S.unset_t0 (); % unset the reference time The input arguments are: WWidth. Particles hit the screen if −W/2≤x≤W/2. [mm] HHeight. Particles hit the screen if −H/2≤y≤H/2. [mm] TTime window. Particles are stored if −T/2≤(t−t0)≤T/2. [mm/c] t0 the reference time [default unset] By default, a new Screen has infinite extension and a boundless time window. 5. Collective Effects In RF-Track, multiple single-particle and collective effects can be added to each element and overlapped. If, for example, one needs to simulate two collective effects, say EFFECT1 and EFFECT2, one can do: E = Drift (1.0); % 1m long drift, or whichever other element E.add_collective_effect (EFFECT1); E.add_collective_effect (EFFECT2); to have both effects act on the beam when the beam travels through the element E. Collective effects in Lattice In Lattice, the effects are uniformly distributed along the element over a user-defined number of steps set through the method set_cfx_nsteps(): E.set_cfx_nsteps (NUMBER_OF_STEPS); The kicks due to collective effects are applied using a velocity Verlet algorithm, which is similar to the leapfrog method. This means that, for example, if only one integration step is used, the algorithm will apply a half-kick at the entrance of the element, perform the tracking through the entire element in one step, and then apply the second half-kick at the exit. Collective effects in Volume Within a Volume, the user must specify how frequently collective effects should be computed by setting the tracking parameter cfx_dt_mm , which defines a time step expressed in mm/ c . Typically, cfx_dt_mm is larger than the integration step dt_mm used to solve the equations of motion. Here is an example: 102 Chapter 5. Collective Effects V = Volume(); V.dt_mm = 0.1; % mm/c, fine integration step V.cfx_dt_mm = 10; % mm/c, apply collective effects every 10 mm/c Collective effects in RF-Track are computed using parallel algorithms. However, the degree of parallelism achievable when tracking without collective effects is inherently superior, as particles are fully independent and there is no interaction between them. For this reason, simulations run faster when cfx_dt_mm is large compared to dt_mm . Of course, there is a trade-off between simulation speed and the accuracy of the results. A convergence study is therefore always recommended. 5.1 Space charge 103 5.1 Space charge Space charge effects are handled differently from other collective effects. While space charge is inherently associated with the bunch itself, other collective effects typically depend on specific components of the beamline. For instance, wakefields arise from accelerating structures (e.g., shortrange wakes) or from particular sections of the beamline (e.g., beam pipes producing resistive-wall wakes). Likewise, spurious magnetic multipole fields are confined to specific magnetic elements. 5.1.1 Space-charge in Volume As discussed in Chapter 2 (Beam Models), the Volume (time integration) environment is the most appropriate for space-charge-dominated regimes, such as injectors and other low-energy systems. In Volume , space charge effects are modeled as acting continuously throughout the entire 3D space. Space charge can be activated by setting the tracking parameter sc_dt_mm , which specifies the kick interval, i.e., how often a space-charge kick is applied during tracking. The unit is mm/c: V = Volume() V.sc_dt_mm = 10; # mm/c As with cfx_dt_mm , it is generally advantageous to choose a space-charge kick interval sc_dt_mm that is larger than the integration step dt_mm . This improves computational efficiency by enabling highly parallel tracking, particularly in situations where the bunch charge distribution remains nearly constant between successive kicks. 5.1.2 Space-charge in Lattice Space charge effects can also be simulated in the Lattice (space integration) environment, but additional considerations are required due to the nature of space-based tracking. When using space integration, i.e., with a Bunch6d , all particles are located on the same longitudinal plane at each step, and their arrival time is represented by the fifth phase space coordinate. However, the evaluation of the space charge force requires the spatial configuration of particles at the same time, rather than their time of arrival at a common longitudinal location. To address this, RF-Track performs an on-the-fly transformation of the bunch: the particles are virtually drifted so that they become simultaneous in time, with their average longitudinal position aligned to the current integration position. This transformation is used solely for the computation of the space charge kicks. A drift step is applied to each particle to achieve this time alignment. Once the space charge kicks are computed and applied, tracking proceeds as usual. This technique ensures that space charge effects can be accurately modeled within the Lattice framework, despite its inherently spatial nature. In the Lattice environment, space charge effects are enabled by specifying the number of space charge kicks to be applied within each element. This is done on an element-by-element basis: E = Element(); % any element E.set_sc_nsteps(NUMBER_OF_STEPS); 104 Chapter 5. Collective Effects Here, NUMBER_OF_STEPS is a positive integer indicating how many evenly spaced space charge kicks should be applied across the element’s length. A larger number of steps increases accuracy but may increase computation time. 5.1.3 Space Charge Models RF-Track provides two algorithms for modeling space charge effects: SC = SpaceCharge_PIC_FreeSpace(Nx=16, Ny=16, Nz=16); SC = SpaceCharge_P2P(); SpaceCharge_PIC_FreeSpace() implements a free-space fast Fourier solver based on 3D retarded integrated Green’s functions. This is the default and recommended space charge algorithm. The parameters Nx , Ny , and Nz define the number of mesh points used along each spatial dimension. These values can be any positive integers—odd or even—and do not need to be powers of two. The space charge field is recomputed on the updated grid at every kick. SpaceCharge_P2P() implements a particle-to-particle algorithm that computes the electromagnetic interaction between all pairs of macroparticles. This method is significantly slower and more susceptible to unphysical artifacts, especially when particles with large charge are too close to one another. It was originally implemented for academic purposes and retained for consistency checks. Use of SpaceCharge_PIC_FreeSpace() is strongly recommended in all practical scenarios. Importantly, unlike many other simulation codes, RF-Track accounts for both electric and magnetic contributions to the space charge force. As a result, beam–beam interactions are naturally included in the simulation. Setting the space-charge algorithm The space charge computing engine can be changed by setting the RF-Track settings’ variable SC_engine. In Octave, this is done like in the following example: SC = SpaceCharge_PIC_FreeSpace(40, 40, 40); RF_Track.SC_engine = SC; In Python, SC = SpaceCharge_PIC_FreeSpace(40, 40, 40) RF_Track.cvars.SC_engine = SC This is usually done at the beginning of a script, as it affects all subsequent space charge calculations. The default engine is SpaceCharge_PIC_FreeSpace(32, 32, 32). Space-charge force The space charge routines can also compute the space-charge force outside of tracking using the method compute_force(). For example, B0 = Bunch6d ( ... ); % A bunch, or a Bunch6dT SC = SpaceCharge_PIC_FreeSpace (40, 40, 40); % SC 5.1 Space charge 105 F = SC.compute_force (B0); % returns the force in MeV/m F is a 3-column matrix, with as many rows as particles in the bunch B0 . The three columns are Fx , Fy, and Fz, the three components of the force exerted on each particle, expressed in MeV/m. 112 Chapter 5. Collective Effects 5.6 Wakefields RF-Track implements three wakefield models: ShortRangeWakefield , which implements K. Bane’s approximation described in [SLAC-PUB-9663, 2003], LondRangeWakefield for multibunch simulations, and Wakefield_1d for user-defines shortand long-range wakefield function provided by the user. 5.6.1 Short-range wakefield ShortRangeWakefield implements the transverse and longitudinal effect of short-range wakefields. SRWF = ShortRangeWakefield (a, g, l); SRWF = ShortRangeWakefield ([ai af], [gi gf], [li lf], Ncells); In its second form, this effect can simulate cell-to-cell variations such as the tapering of iris apertures. The required input parameters describe the accelerating structure’s cell geometry: athe average iris aperture radius [m] gthe average gap length [m] lthe average cell length [m] [ ai af ] a 2-element vector with the initial and final iris apertures [m] [ gi gf ] a 2-element vector with the initial and final gap lengths [m] [ li lf ] a 2-element vector with the initial and final cell lengths [m] Ncells The number of cells in the structure [INTEGER] The meaning of these quantities is shown in Figure 5.1. In evaluating the wakefield effect, W⊥(s) and W||(s) are calculated using the analytic approximation of the normalized wake potential presented in [K.Bane, SLAC-PUB-9663]: w⊥(s) = 4Z0cs⊥0 πa41−1+rs s⊥0exp−rs s⊥0 [V/pC/m/mm] wk(s) = Z0c πa2exp−rs sk0[V/pC/m] where Z0 is the impedance of free space, a is the average aperture radius of the structure, g is the gap length, dis the length of the cell, sk0and s⊥0are sk0=0.41 a1.8g1.6 d2.5, s⊥0=1.69 a1.79g0.38 d1.17 , respectively. The range of validity of these formulægenerally covers most of the cases. Inquiring the wake function The user can inquire RF-Track about the single-particle Wake function using the method: Wl = SRWF.w_long ( S ); % returns V/pC/m Wt = SRWF.w_transv ( S ); % returns V/pC/m/mm where the input parameter: S is the longitudinal distance in meters from the source particle. Since the Wakefield can only affect the trailing particles, S<0. [m] 5.6 Wakefields 113 g l a Figure 5.1: The geometric parameters a,g, and lused to describe the short-range wakefield. Cell-to-cell misalignment To study the impact of cell-to-cell misalignment, RF-Track can randomly scatter each cell using the method: SRWF.scatter_cells (RANGEX, RANGEY ); % mm The parameters RANGEX and RANGEY determine the amplitude of the offset in mm. Each call to this function offsets each of the Ncells cells of the structure within a range −RANGE/2 and +RANGE/2 according to a uniform random distribution. 114 Chapter 5. Collective Effects 5.6.2 Long-range wakefield In the case of multi-bunch operation, long-range wakefield effects can become significant. In RF-Track, the command LongRangeWakefield allows the simulation of long-range wakefield effects. Its input arguments are: LRWF = LongRangeWakefield (A, freq, Q, angle); LRWF = LongRangeWakefield (A, freq, Q); The required input parameters describe the transverse high-order modes: AA vector or a matrix containing the modes’ amplitudes [V/pC/mm/m] freq A vector or a matrix containing the modes’ frequencies [GHz] QA vector or a matrix containing the modes’ Qfactors [unit less] angle A vector or a matrix containing the modes’ polarization angles. If NaN , no polarization is considered, and the kick is given in the xand ydirections. (DEFAULT: NaN) [deg] When A,freq, and Qare matrices, e.g. A=   A11 A12 ··· A1m . . .. . . An1An2··· Anm   ,(5.1) the first index, n , runs over the longitudinal axis of the structure, whereas the second index, m , runs over the different modes. Thus, one can simulate a set of modes that varies along the structure. The applied normalized wakefield is: w⊥(s) = −∑ n An·expπs Qnλnsin2πs λn,[V/pC/m/mm] wk(s) = 1 2π∑ n λnAnexpπs Qnλncos2πs λn,[V/pC/m] with λn being the wavelength corresponding to the frequency fn . The polarization angle indicates the polarization of the mode. If angle=0 , for example, the wakefield will only affect the x axis. If angle=90 , for example, the wakefield will only affect the y axis. If angle=NaN , the wakefield will equally affect the xand yaxes. Referring to the notation in literature, one can write the coefficients Anas follows: An=cRs/Q 2a2 1 Lstr =cZ/Q 2a2,[V/C/m/mm] with Rs the shunt impedance in Ohm, a the source particle’s offset, Lstr the structure length, and Z the normalized shunt impedance in Ohm/m: Z=Rs Lstr .[Ω/m] (See Wangler). Inquiring the wake function The user can inquire RF-Track about the single-particle Wake function using the method: 5.6 Wakefields 115 Wl = LRWF.w_long ( S ); % returns V/pC/m Wx = LRWF.w_x(S);% returns V/pC/m/mm Wy = LRWF.w_y(S);% returns V/pC/m/mm where the input parameter: S is the longitudinal distance in meters from the source particle. Since the Wakefield can only affect the trailing particles, S<0. [m] If no polarization is specified, w_x(s) and w_y(s) return the same value. 116 Chapter 5. Collective Effects 5.6.3 Generic wakefield The collective effect Wakefield_1d allows the user to provide an arbitrary single-particle wake function, or Green’s function, using two 1d vectors sampling the function. WF = Wakefield_1d (Wt, Wl, hz); The required input parameters to describe the wakefield are: Wt A 1d vector with the transverse component of the wake function [V/pC/mm/m] Wl A 1d vector with the longitudinal component of the wake function [V/pC/m] hz the 1d mesh spacing of Wt and Wl [m] Inquiring the wake function The user can inquire RF-Track about the single-particle Wake function using the method: Wl = WF.w_long ( S ); % returns V/pC/m Wt = WF.w_transv ( S ); % returns V/pC/m/mm where the input parameter: S is the longitudinal distance in meters from the source particle. Since the Wakefield can only affect the trailing particles, S<0. [m] 5.7 Passage of particles through matter 117 5.7 Passage of particles through matter Three complementary effects can be turned on or off independently to simulate the passage of charged particles through matter. The first is multiple Coulomb scattering; the second is stopping power, which accounts for the energy loss described by the Bethe-Bloch formula; and finally, energy straggling, which is the fact that the energy loss isn’t the same for all particles. 5.7.1 Multiple Coulomb scattering The Multiple Coulomb Scattering module allows the simulation of the interaction between the bunch particles and a material. This effect can be added to any elements to simulate beam tracking through matter. Follows a list of the constructor and the main methods: M = MultipleCoulombScattering( material ); M = MultipleCoulombScattering( X0, Z, A, density, I=-1 ); M.enable_log_term(); M.enable_fruehwirth_model(); % DEFAULT M.enable_wentzel_model(); % DEFAULT M.disable_log_term(); % DEFAULT M.disable_fruehwirth_model(); M.disable_wentzel_model(); The input arguments describing the material are: material a pre-defined material among: ‘ air ’, ‘ water ’, ‘ beryllium ’, ‘ lithium ’, ‘ liquid_hydrogen ’, ’titanium’, ‘tungsten’ [STRING] X0 the radiation length [cm] Zthe atomic number of absorber [INTEGER] Athe atomic mass of absorber [g mol−1] density the material density [g cm−3] I the mean excitation energy. If not provided, RF-Track will compute it using the approximation proposed by Bloch. [eV] 5.7.2 Stopping power Charged particles lose energy as they pass through matter, governed by the Bethe-Bloch equation, which describes how the energy loss depends on the properties of the material and the particle. The energy loss per unit path length is defined as −dE ds =4πNA α2 me ρZ A 1 β2h1 2ln 2meγ2β2Tmax I2−β2i.(5.2) In RF-Track, Eq. (5.2) is automatically applied as the default energy loss model when the command MultipleCoulombScattering(material) is used. For electrons in water and electrons in air, RF-Track uses the data tabulated by the National Institute of Standards and Technology, https://physics.nist.gov/PhysRefData/Star/Text/ESTAR.html , for muons in liquid 118 Chapter 5. Collective Effects hydrogen, RF-Track uses the data tabulated by the Particle Data Group https://pdg.lbl.gov/ 2022/AtomicNuclearProperties/MUE/muE_hydrogen_liquid.txt. 5.7.3 Energy straggling While a bunch traverses a material, the energy loss is not the same for all the particles. Hence, particles incident on the foil with the same energy emerge with different energies. This energy distribution caused by the scattering material is known as the particles’ ’straggling”. 6. Inverse Compton Scattering 6.1 Introduction RF-Track can simulate laser beams. The laser can interact with electrons, positrons, or any other charged particle through Thompson and Compton scattering. The element LaserBeam allows the laser beam’s initialisation. RF-Track provides several functions that can fine-tune the laser-beam interaction conditions. Position and angle offsets can be applied to the laser beam, enabling the study of misalignments and other imperfections. 6.1.1 General parameters Structure LaserBeam contains all information required to define a laser beam from its general parameters. The possible arguments are: LB = LaserBeam(); LB.pulse_energy; % mJ, laser pulse energy LB.pulse_length; % ps, laser pulse length LB.wavelength; % nm, laser wavelength LB.length;% m, length of the interaction region LB.P; % laser polarization, abs(P)<=1 or NaN for unpolarized laser beams LB.R; % mm, laser sigma spot size, can be Gaussian or tophat profile LB.Rx; % mm, horizontal sigma spot size LB.Ry; % mm, vertical sigma spot size LB.M2 = 1; % laser beam quality factor LB.rep_frequency = 0; % Hz, the laser repetition frequency LB.number_of_pulses = 1; % number of pulses in the train The laser beam size can be defined in two ways: for a symmetric laser beam profile, you can use 120 Chapter 6. Inverse Compton Scattering LB.R ; if the beam is asymmetric along two orthogonal axes, you can use LB.Rx and LB.Ry to define the two transverse sizes. Note that LB.R is half the 1/e2laser waist radius, w0. The parameter LB.P allows you to set the polarisation of the laser beam. At the moment, only linear polarisation has been implemented. If the unit vector ˆ k defines the polarisation axis of the incoming laser beam in the plane orthogonal to the direction of propagation of the laser, then ˆ k=  cosθ sinθ 0 , and P=sinθ. The default setting, LB.P = NaN, means unpolarised laser beam. Note that the length of the interaction region LB.length can be zero, i.e., an interaction point can be inserted as a thin element in a lattice. Even in the case of a zero-length interaction region, RF-Track calculates the laser-beam interaction by reconstructing the full 3D shape of both bunch and laser pulses from the moment of first contact until the moment of separation of the two colliding bunches. Define the direction of the laser The direction of propagation of the laser beam is completely arbitrary and is specified as a 3D vector. LB.set_direction(nx, ny, nz); The possible arguments are: nx horizontal component [a.u.] ny vertical component [a.u.] nz longitudinal component [a.u.] For example, for a full head-on collision, one must set nx=0 , ny=0 , and nz=−1 . Notice that the vector defined by nz , ny , and nz doesn’t need a unitary length; whichever length, only its direction is considered. 6.1.2 Particles-laser interaction Distinct functions and parameters define the laser-beam interactions in the LaserBeam element. Define the interaction point If the length of the LaserBeam element is not zero, you need to specify the position of the interaction point between the incoming beam and the laser along the longitudinal axis. LB.set_position(Z); The parameter Zmust be a number between 0 and the length of the element. Set number of generated macro-particles per slice RF-Track implements macro-particles for tracking. The beam(s) generated from interactions are also structured in macro-particles. This function allows you to adjust the number of macro photons generated in each slice to ensure enough macro photons are created in one run to generate statistics. This affects both the precision of the interaction simulation and its runtime. 6.2 Collecting the generated photons 121 LB.min_number_of_gammas_per_slice(nr_gamma); Where nr_gamma is an integer for the number of photon macro-particles simulated per slice. Define the number of steps taken by the simulation Tracking through elements is inherently discrete, as the beam is advanced through an element in a discrete number of steps. The number of steps taken across the interaction space can be changed in RF-Track. A large number of steps improves the precision of the simulation but increases the runtime. LB.set_nsteps(nsteps); Where nsteps is an integer, the number of steps taken during the simulation across the interaction region. An odd number of steps is recommended to capture the interaction in the central step. 6.2 Collecting the generated photons When a bunch traverses a LaserBeam element, RF-Track computes the Compton scattering between the bunch and the laser. All photons generated during the interaction are added to the travelling bunch. The user can retrieve them using the Bunch6d ’s method get_phase_space() . In particular, inquiring about the mass or the charge of the particles in the bunch allows the discrimination between the bunch’s particles and photons. Follows an example: % L is a lattice containing a LaserBeam element B1 = L.track(B0); % Read the output phase space M1 = B1.get_phase_space(’%x %xp %y %yp %t %P %m %N’); % Separate the photons from the electrons is_photon = M1(:,7) == 0; % pick the photons (mass == 0) M1p = M1( is_photon,:); % photons M1e = M1(~is_photon,:); % other particles % Total number of single photons generated Np = sum(M1p(:,8)); % Plot the XX’ photons’ phase space scatter(M1p(:,1), M1p(:,2)); xlabel(’x [mm]’); ylabel(’x’’ [mrad]’); 128 Chapter 8. Extending RF-Track 8.2.2 Implementation To define a UserElement , one creates a subclass of RF_Track.UserElement and overrides the appropriate method to apply the desired transformation to the particle coordinates. 8.2.3 Example: A Custom Element in Python The following example defines a simple user element that transports the bunch (in this case, a single particle) and appends a new particle to it. All user-defined elements must inherit from the Python base class RF_Track.UserElement and implement a method: UserElement.track(self, bunch) The required input arguments are: self The instance of the user-defined element class. bunch The bunch being tracked through the element (type: Bunch6d ). The custom Python class must initialize the parent class by calling its __init__() method, passing the element length (in meters) as argument. The input bunch can be freely modified within the track() method. New particles can be added using the method bunch.append(new_bunch) , where new_bunch is an instance of Bunch6d or an extended phase-space, as described when illustrating the Bunch6d constructors, in section 2.2.1. import RF_Track as rft import numpy as np class myElement(rft.UserElement): def __init__(self, length=0.0): super().__init__(length)# Passing the length [m] is mandatory def track(self, bunch): print(’PYTHON::Hi from Python!’) print(f’PYTHON::bunch.S = {bunch.S} m\n’) # Modify the phase space of the incoming bunch bunch[0].x += 5 bunch[0].Pc = 50 # New particles are created through a new phase space # np.array([ %x %xp %y %yp %t %P %mass %Q %N ]) new_particles = np.array([ 1, 2, 3, 4, 5, 6, rft.muonmass, +1, 1e4]); # Add new particles to the current bunch 8.3 Custom Electromagnetic Fields with UserField 129 bunch.append(new_particles) # User-defined element P = myElement(length=1.0) # Crate a lattice containing it L = rft.Lattice() L.append(P) # Create one test particle and track it phase_space = np.array([ 0., 0., 0., 0., 0., 100. ]) B0 = rft.Bunch6d(rft.electronmass, 0.0, -1, phase_space) B1 = L.track(B0) # Print results (two particles) print(B1.S) print(B1.get_phase_space(’%m %x %xp %y %yp’)) # mass x xp 8.3 Custom Electromagnetic Fields with UserField The UserField construct enables the definition of user-defined electromagnetic fields—either static or time-dependent—which can be attached to a Volume or a Lattice element. 8.3.1 Applications • Implementing analytical field models. • Interpolating custom field maps from measurements or simulations. • Simulating pulsed magnets or time-varying RF fields. 8.3.2 Field Evaluation To insert a user-defined field, the user needs to define a method to return the electromagnetic field at a given point in space and time. This method is called internally by RF-Track when integrating the equations of motion, or when using Volume::get_field() for plots. Fields must have a finite length larger than zero in order to be properly inserted in a RF-Track’s Volume or Lattice. The length must be passed to the parent class through its method __init__(length). Note: Due to the complexity of integrating a user-defined field into the core tracking engine of RF-Track, simulations involving a UserField must run in single-threaded mode. As a result, RF-Track cannot utilize multiple CPU cores in this case. To enforce this, the number of threads must be explicitly set to one at the beginning of each script that uses a UserField: import RF_Track as rft rft.cvar.number_of_threads = 1 130 Chapter 8. Extending RF-Track This limitation may significantly slow down complex or large-scale simulations. Users are advised to isolate user-defined fields to the minimum necessary sections of the beamline. 8.3.3 Example: A Custom Electromagnetic field in Python The following snippet defines a simple collective effect that applies a constant transverse kick to all particles. import RF_Track as rft import numpy as np rft.cvar.number_of_threads = 1 # see text # Define a user-supplied quadrupole magnetic field class myQuadrupole(rft.UserField): def __init__(self, length=0.0, gradient=0.0): super().__init__(length)# Passing the length [m] is mandatory self.G = gradient / 1e3 # Convert from T/m to T/mm def get_field(self, x, y, z, t): print("PYTHON::Hi from Python!") # Magnetic field for ideal quadrupole: Bx = G*y, By = G*x Bx = self.G *y By = self.G *x Bz = 0.0 # No electric field E = np.zeros(3) B = np.array([Bx, By, Bz]) return E, B # Instantiate the custom quadrupole field with gradient 4 T/m Q = myQuadrupole(length=0.2, gradient=4.0) # Add the field to a Volume, at position (0, 0, 0) V = rft.Volume() V.add(Q, 0.0, 0.0, 0.0) # Query the field at position (x=0, y=2, z=4) mm and time t=6 mm/c E, B = V.get_field(0.0, 2.0, 4.0, 6.0) # Print the field values print(f"E field = {E} V/m") print(f"B field = {B} T") 8.4 Custom collective effects with UserEffect 131 8.4 Custom collective effects with UserEffect The UserEffect construct allows users to define custom collective effects that act on the particle distribution. These effects can be attached to any RF-Track element and evaluated at configurable intervals along the element length. 8.4.1 Capabilities • Applying custom forces based on the beam distribution. • Modeling beam–environment interactions (e.g., wakefields). • Introducing energy loss mechanisms or external feedbacks. 8.4.2 Returning modified bunches The user-provided function may optionally return a modified bunch. This enables the implementation of complex effects such as energy spread growth, halo generation, or particle-matter interactions with secondary production. 8.5 Secondary particle generation and external interfaces The three constructs — UserElement , UserField , and UserEffect —, used at runtime during tracking, support the injection of new particles into the simulation. This feature is essential for coupling RF-Track with external codes, for example: •Geant4 — for detailed simulation of particle–matter interactions. •BDSIM — for radiation and shielding studies. •FLUKA — for hadronic cascade modeling and energy deposition. New particles can be added by returning a new Bunch6d or Bunch6dT object from the userdefined function. These secondaries are then added to the bunch and tracked like all other particles, interacting with fields and collective effects. 8.5.1 Example: A Custom Collective Effect in Python All user-defined collective effects must inherit from the Python base class RF_Track.UserEffect and implement the method compute_force(). [force, new_particles] = UserEffect.compute_force(self, bunch, S_mm, dS_mm) The required input arguments are: self The instance of the user-defined collective effect class. bunch The bunch being subjected to the force. Type: Bunch6d or Bunch6dT. S_mm The initial position along the element [mm]. dS_mm The step length over which the force acts on the bunch [mm or mm/cin case of Bunch6dT]. The method must return two values: 132 Chapter 8. Extending RF-Track force An N×3 array, where N is the number of particles in bunch . Each row contains the components of the applied force (Fx,Fy,Fz)in units of MeV/m. new_particles A bunch containing any new particles to be added to the simulation. Type: Bunch6d or Bunch6dT. Can be empty or omitted. Note: The incoming bunch must not be modified directly. Collective effects are intended as thin kicks. Although the parameter dS_mm must always be provided (and may be used internally by the effect—for example, when computing scattering events or radiation emission). Any new particles should be returned through the new_particles object instead. The following code defines a simple collective effect that applies a constant transverse kick to all particles: import RF_Track as rft import numpy as np class mySpaceCharge(rft.UserEffect): def compute_force(self, bunch, S_mm, dS_mm): # Apply a uniform transverse kick N = bunch.size() # number of particles in the bunch # Compute the force acting on the existing particles force = np.zeros((N, 3)) # Fx, Fy, Fz force[:,0] = 0.1 # Kick in x [MeV/m] force[:,1] = 0.1 # Kick in y [MeV/m] # Add new particles new_particles = rft.Bunch6d([...]) return force, new_particles This minimal implementation can be extended to include diagnostics, conditionals based on beam shape, or coupling with external codes or data. The optional second return value (here, the unchanged bunch) may be used to append new particles to the beam. 8.6 Traversing Beamlines with UserVisitor RF-Track provides a user-extensible visitor interface that allows custom traversal of a Lattice or a Volume . This is implemented via the UserVisitor class, which follows the classical visitor pattern. It enables Python users to define actions that are applied to each element of a beamline during its traversal. 8.6.1 Key use cases UserVisitor can be used to perform read-only inspections, logging, graphical rendering, or even to build auxiliary data structures. Some example applications include: • Generating a survey report of all elements and their 3D coordinates. 8.6 Traversing Beamlines with UserVisitor 133 • Drawing a schematic of the beam line for visualization. • Exporting lattice data to an external format (e.g., JSON, XML). 8.6.2 Implementing a Visitor in Python To define a custom visitor, the user must subclass RF_Track.UserVisitor and implement the method: def visit(self, element): ... This method is called by RF-Track during traversal. The parameters are: element The current item being visited. Type: Element or any derived classes. The user can detect the exact type of each element, using the following construct: if type(element) == rft.Quadrupole: print(’The element is a Quadrupole’) or alternatively if isinstance(element, rft.Quadrupole): print(’The element is a Quadrupole’) A visitor can access and process each object’s parameters, fields, or identifiers. The visitor can stop the traversal early by calling the method: my_visitor = myCustomVisitor() my_visitor.stop_at(element) where element is the desired end element in the lattice. Using the Visitor Once defined, the user visitor is applied to a Lattice or a Volume using the accept() method. For example: L.accept(my_visitor) RF-Track will traverse the internal structure of the lattice (or volume) and invoke the appropriate method on the visitor for each encountered object. Example: Printing Element Locations The following Python example prints the global position of each element in a lattice: 134 Chapter 8. Extending RF-Track import RF_Track as rft import numpy as np ## Define a FODO cell ... ## Start UserVisitor class myVisitor(rft.UserVisitor): def visit(self, e): print(’PYTHON::Hi from Python!’) if type(e) == rft.Quadrupole: g = e.get_gradient() print(f’PYTHON::This is a quadrupole with gradient = {g} T/m’) if type(e) == rft.Lattice: n = e.size() print(f’PYTHON::This is a lattice containing {n} elements’) V = myVisitor() FODO.accept(V) This mechanism provides a clean and powerful way to extend RF-Track’s beamline traversal logic without modifying the core C++ code or duplicating structural logic. 8.7 Getting started Users are encouraged to begin with the examples provided in the RF-Track Jupyter tutorials directory, and consult this manual’s reference sections for further customization details. 9. Examples 9.1 Example of bunch creation 9.1.1 Bunch6d from an arbitrary user-defined distribution For an example of a bunch created from the Twiss parameters, see Chapter 1. The following example creates a bunch from a matrix. The matrix has dimensions N×6 , where 6 is the size of the phase space, and Nis the number of macroparticles in the bunch. % create a bunch of 1e12 100 MeV/c electrons, using 1000 macroparticles P = 100; % MeV/c total momentum O = zeros(1000,1); % define a column vector of zeros I = ones(1000,1); % define a column vector of ones X = randn(1000,1); % mm, column vector of Gaussian-distributed positions Y = randn(1000,1); % mm, column vector of Gaussian-distributed positions % create the beam matrix M=[XOYOOP*I];% Bunch6d %x %xp %y %yp %t %P % create a bunch B0 = Bunch6d(RF_Track.electronmass, 1e12, -1, M); % retrieve the phase space following MAD-X convention T0 = B0.get_phase_space("%x %px %y %py %Z %pt"); % retrieve the phase space following the TRANSPORT convention T0 = B0.get_phase_space("%x %xp %y %yp %dt %d"); 136 Chapter 9. Examples % retrieve the phase space following PLACET convention T0 = B0.get_phase_space("%E %x %y %dt %xp %yp"); % save on disk B0.save(’my_bunch.rft’); % save the beam in RF-Track binary format B0.save_as_dst_file(’my_bunch.dst’, 750.0); % save as DST, 750 MHz RF B0.save_as_sdds_file(’my_bunch.sdds’,’my useful comment’); % save as SDDS % save as an Octave matrix save -text my_bunch.txt T0;  9.1.2 Chirped Bunch6d from Twiss parameters 1%% Load RF-Track 2RF_Track; 3 4%% Beam 5Pref = 100; % MeV/c, reference momentum 6Q = -1; % electrons 7Nparticles = 10000; % 10k particles 8 9%% Twiss parameters 10 Twiss = Bunch6d_twiss(); 11 Twiss.emitt_x = 1; % mm.mrad normalised emittance 12 Twiss.emitt_y = 1; % mm.mrad 13 Twiss.beta_x = 1; % m 14 Twiss.beta_y = 1; % m 15 Twiss.sigma_t = 1; % mm/c, bunch length 16 Twiss.sigma_pt = 5; % permille, energy spread 17 Twiss.disp_z = 1; % m, longitudinal dispersion (chirp) 18 19 %% Create bunch 20 B0 = Bunch6d(RF_Track.electronmass, 0.0, -1, Pref, Twiss, Nparticles); Index Symbols Bunch6dT .............................25 Bunch6d ...............................21 Lattice ...............................41 Volume ................................41 autophase() ..........................52 get_field() ..........................56 A Automatic matching . . . . . . . . . . . . . . . . . . . . . 54 Automatic synchronization . . . . . . . . . . . . . . . 52 B Backtracking...........................54 Beamloading.........................108 Bunchpersistency.......................35 C Citation................................20 Coastingbeams.........................28 E Elements Absorber, Absorber ..............................81 Adiabatic Matching Device, AdiabaticMatchingDevice ...........83 Beam Position Monitor, Bpm....................................98 Coil, Coil ..................................67 Corrector, Corrector ............................66 Custom Effects, UserEffect ..........................131 Custom Elements, UserElement .........................127 Custom Fields, UserField ...........................129 Drift, Drift .................................60 Electron cooler, ElectronCooler ......................82 LaserBeam, LaserBeam ............................86 Multipole, Multipole ............................79 Pillbox cavity, Pillbox_Cavity ......................78 Quadrupole, Quadrupole ...........................61 Rectangular Bend, RBend .................................65