Transcription
Fermi labESME:TM-12742041.000LONGITUDINAL PHASE SPACE PARTICLE TRACKING - PROGRAM DOCUMENTATIONJ. A. MacLachlanMay 1984
21.TM-1274INTRODUCTIONThe program ESME has been developed to model various manipulations of the beamin a circular accelerator that are effected by controlling the radiofrequency accelerating potentials like, for example, stacking, bunch shaping, bunch rotation, changeof harmonic number, and so on.The advantages of this program over earlier versionsand other programs of which theautho is aware are precise treatment of the non-linearity of the equations of ·motion, a convenient and extensive set of routines forgenerating initial distributions, flexibility in the specification of variation ofparam terswith respect to time or properties of the distribution, and rather completegraphical output.Because it was found rather early in the development that theprogram would continue to evolve, some thought has gone into organizing it so thatsuch additions do not usually require pervasive modifications.This programmingeffort has not, of course, entirely escaped the typical difficulties associated withgenerality and elaboration, viz., proliferation of input parameters and errorshiding in infrequently used sections of the code or freshly introduced in a recentmodification.To confront the first difficulty an effort has been made to requirefew parameters for a simple problem by establishing reasonable defaults and byseparating various program features so that only the relevant ones need be considered.The second difficulty prompts.the following disclaimer:the program is the resultof the author's particular needs and interests so that certain features used infrequently or only in early stages of the development may not have been properlymodified as the program evolved.However, the program structure is sufficientlymodular that such bugs as remain are likely isolated to special functions where theycan be easily found when their effects have been characterized.In the detaileddescription which follows a few features have been noted as being less fully tested.Such dubious elements could be excised from the program with little loss from whathave been claimed as its principal virtues.This writeup has been produced to helpnew users at CERN to get started with the code and is being distributed at Fermilabpartly as imf ormation for possible users and partly to document the basis for manydesign calculations for the TeV I Project.2.THE COORDINATE SYSTEMThe phase space coordinates used in ESME are the azimuth of a particle at thetime that the accelerating potential passes through the synchronous phase and itstotal energy in MeV.The values of the azimuth-TI 8 TIare positive upstream
3TM-1274of what is taken as a single accelerating gap.For those very specialcases for which it is signficant to distribute the RF sources, the program can beused if the segments between RF stations are the same and the stations themselvesare identical.Otherwise tedious program modification would be necessary.A resultof this choice of coordinate system is that a well behaved bunch will generallybe centered about 8 0 [mod h].If there are h equivalent bunches in a calculation,one will generally choose to look at a single one by setting the periodicityparameter (called FRAC in the program) equal to the harmonic number h so that thecalculation is restricted to n/h with periodic boundaries, i.e., what leaves oneboundary immediately reappears at the other.The definition of the phasesw1and 2 of the RF sources can be subtle in somecomplicated prodcedures, but by reference to typical cases the desired phase programsfor the complicated cases can be constructed.The simplest example is a stationarybucket produced by a single frequency source.In this casethe lattice is below or above transition, i.e.,added to change the bucket shape,w1wl w1TIaccording as .If another harmonic h 2 issis unchanged and 2 const. is whateverrelative phase is required to produce the desired bucket shape.some value other than 0 orwlcan be set to , for example, to stretch a bunch along the spearatrix atsthe unstable fixed point.When p# 0, 1 usually is scalculated from vl, but if v2 is present theremay be corrections resulting from acceleration produced by V2 Both land 2maybe changed linearly, individually or in concert, during the initiation or conclusionof a stacking operation, for example.to 2may also change continuously with respectif RF source 2 is being run asynchronously as specified by a momentum error1parameter.Because of the variety of possibilities it has proven helpful to fix RF source1 as the principal one and regard source 2 as in some way auxiliary so that thephase of source 2 is always defined relative to the phase of source 1.To maintainthe primacy of source 1 in a case where one source starts out dominant but becomesless important during the manipulation (e.g. rebunching at a new harmonic), theparameters associated with the sources are interchanged when the bucket heightproduced by source 2 exceeds that produced by source 1.height from source 2 is greater, it becomes source 1.That is, when the bucket
43.TM-·1274THE DIFFERENCE EQUATIONSIn the conventional derivation of the differential equation for synchrotronoscillations one starts from equations like: E.1,n eV sin J .n1,n Pi,n l l i,n wnT i,n lwhere E.is .the energy of the i th particle on its nth passage around the ring,1,n P.is phase of the radiofrequency of angular frequency w and amplitude V when1,nnnith particle makes its nth traversal of the accelerating gap, and l is the1,ncirculation period for the ith particle on its n 1st circuit of the ring. By settingT.a synchronoism relation between the circulation frequency and the radiofrequency andtakinge.1. v./c1 esfor all particles, one arrives at a familiar second ordernon-linear differential equation.This equation is sometimes solved numerically incalculations of behaviour of particle distributions.When recourse is made to acomputer, however, one can achieve a far more accurate treatment of the abovedifference equations by solving them directly, because they are appropriate to thecomputer nearly as they stand.Then one can precisely calculate each Si and therebytreat the kinematics precisely at nearly the same rate of computation as one caniterate the difference equation approximation to the already approximate differential on E . can be included throughn,in,1any order to which the dependence of the lattice parameter yt (transition energy gamma)equation.The nonlinearity in the dependence ofTon momentum is known.The treatment starts from the assumption of a reference orbit of known averageradius Ron which the average field B is known.yreference orbit has a momentum:np0 2.997925A particle which follows thex 10 2 B Ry0where p 0 is in MeV/c, Ro is in meters and B is in Tesla.yPo hasso that the angular frequency of its circulation isA particle of momentum
TM-12745If one synchronizes the RF frequency on each turn so thatw hrlno,n- hrls,nthen Pi,n lcp,1. ,n 2Tih rls,nQi,n lTo simplify calculations in which there is more than one RF frequency present it isconvenient to replace the RF phase variable by an azimuthal angle 8 measured from theRF source location in a direction opposite to the direction of particle circulation.e(This convention disposes of the usual - sign relatingto P);8, P./h.i.,ni.,nIfn.1. ,n l Qs ,n lturn.then the cj -equation shows that 8 as defined will increase by 2TI eachThus it is natural to write the cj -equation in a "rotating frame" by subtracting 2TIeach turn.With these changes the difference equations become: E.i.,n eVA 1 )sin(h 1 8. i.,nq/n1» eVn( 2 )sin(h 2 8i,n .Qe.i.,n l e.1,n 2TIs,ni,n 1Q,- 1where the possibility of multiple radiofrequencies each with its own amplitude andindependent phase has been indicated.In fact, because the rhs of the energy equationis simply a Fourier expansion of an arbitrary periodic function it is obvious that anytype of repetitive waveform can be used; the sinusoidal form is merely a familiar caseresulting from the natural properties of resonant RF sources.As disucssed in thepreceding section on coordinates, the definition of phase used in the program has been(8, E) 1 cjis' the synchronous particle is at(O,Es)' and 2 is defined relative to the zero of sin(h 1 8 1 ).specialized slightly so that in a typical case
64.TM-1274DATA4.1OverviewGenerally input data are read according to NAMELIST statements so thatonly parameters which need to be changed from current values need to be provided.Furthermore, because defaults have beeiovided for many parameters which areappropriate in typical applications, a simple problem can be specified with verylittle data.Although the defaults do not provide a complete input set, there aredefaults even for parameters which are usually specific to a given problem.Thesedefaults have been provided so that neglect to specify a parameter will not necessarily cause the program to fail; if the program continues with a default value inplace of the appropriate value the results will, of course, be wrong, but the factthat execution continues may permit more than a single error to be found in a testrun.The input parameters are nearly all distributed by function amongst sevenNAMELISTs in subroutines serving those functions.The subroutines are called inresponse to single character commands in the data set which invoke the desired programfacilities.These data commands are discussed in the first section following, andthe contents of the various NAMELISTs are detialed in succeeding sections.The onlydatum routinely used which is not NAMELIST input is a title for the graphical output.The program has a number of dummy enteries for the installation of special purposeroutines.Such routines may request additional parameters and these may well not be inthe NAMELIST format as the user prefers.It should be noted that for many paraemtersthe order of entry is immaterial, but it is important that the lattice parameters beentered (R command) before the RF parameters (A command) and that both be enteredbefore establishing the initial particle distribution (P command).Thus a normaldata sequence will be R, A, P with other connnands in any convient order.The order ofdisucssion of the commands reflects this important point.4.2Summary of the Data CommandsThe principal functions of the program are activiated by single character commandsin the data.These commands may lead to the reading in of data relevant to a particularaspect of the problem and/or action upon the data.commands will be detailed n later sections.purposes are tabulated.The data requested by particularIn this section the commands and theirThere is some mnemonic value in the letter associated with agiven function, and the intended association is indicated by underlining a letter inthe command purpose.Underling of the command itself is meant to indicate that thisis a basic command likely to be needed in any application.Remember the R, A, and P
TM-12747are most important and should be invoked in that order.COMMANDFUNCTIONread in lattice (E ing) parameters and initialize quantities dependentsolely on lattice.read RF (acceleration) paraemters.populate phase space with an initial distribution.read in parameters governing the tracking cycle.read in parameters governing graphical output.execute particle tracking according to current paraemters.s uit the program display the distribution (periodic displays will occur without D command byconditions specified with 0).retreive (get) a set of parameters from TAPE7 (if all are good, no need forR, C, A, P-; and 0).!!.ave distribution and current parameters on TAPE8.enter dummy routine SHAZAM at main entry to .E inker with any quantity inCOMMON storage - almost anything except loop counters and arithmeticalin'termediates.additional entries to SHAZAM to provide arbitrary user defined functions.RApc0E.9.DGsT19wwrite a comment into the printed output.initialize the random number generator with an octal integer read accordingto 02 format.set up and execute history plots of selected parameters.calculate the FFT of the 8 distribution and set up for analysis every turnif desiredIHF4.3Lattice ParametersThe parameters are read in NAMELIST \RING\ and the values stored in COMMON \RINGP\.This read is in SUBROUTINE RINGPAR which also derives quantities liken, yt' etc. whichdepend only on lattice parameters.RINGPAR is called by the data command R.Theinput parameters are:NAMEMEANING OR FUNCTIONUNITREQGAMMATThe average radius of the central orbitTransition gamma (ALPHA )- 2ALPHA ICoefficient of3 p/pin a i:P i OALPHA2ALPHA3EKALPHa. ( p/p)miDEFAULTNoneNoneNone1where p/p is taken w/central orbit unless EKALPH # 0is given as kinetic energy of the orbit to which a.1are referred2Coefficient of ( p/p) in series for a3Coefficient of ( p/p)in series for app.Kinetic energy corresponding to the momentum of .theorbit about which a is expandedpMeVEK
TM-12748EK Kinetic Energy of the central orbitMeVNoneEK2Kinetic energy at the end of a magnetic field change(EK2 O. B const.)MeVo.TlStarting time for magnetic field change(Tl 0 start tracking part way through ramp)so.T2End time for magnetic field change operative for EK2sNoneKURVEBSpecifies form of magnetic field ramp1.2.3.4.FRAC101Ilinear ramp from P (EK ) to P2 (EK2)increasing parabolic rampbiased sinusoidal rampdecreasing parabolic rampAzimuthal periodicity; limits calculation to- /FRAC THETA /FRAC with periodic b.c.14.4 Radiofrequency ParametersThe parameters specifying the accelerating potential are read according toNAMELIST jRFj which resides in SUBROUTINE RFPROG which is invoked by the data command A.The input parameters are:NAMEFUNCTION OR MEANINGHlHarmonic number for the principal RF source1.H2Harmonic number for an auxiliary RF source, if any.Program will force the source producing the greaterbucket height to be source 11.HWl·Harmonic number of the "supression" of source 1(for suppressed bucket systems; HWl 1 is normal)1.HW2Harmonic of bucket suppression for source 21.PSilIInitial phase of RF source 1rado.PSI2IInitial phase of RF source 2rado.PSilFFinal value of phase of source 1 when phase changehas been requested (linear variation provided)rado.PSI2FFinal value of phase of source 2radO.MPHlNumber of beam turns for change from PSilI to PSilF1MPH2Number of beam turns for change from PSI2I to PSI2F1Vl Initial peak voltage of RF source 1V2 Initial peak voltage of RF source 2VlFinal peak voltage of RF source 1V2Final peak voltage of RF source 2UNITMVDEFAULT1.xl0- 6
TM-12749MlNumber of beam turns for change from Vl0 to Vl1M2Number of beam turns for change from V20 to V21DELTA2Momentum off set p/p corresponding to the frequencydifference between source 2 and h 2n . Source 1 isalways synchronous. (This is a lit le used feature andmay have problems.)1KURVE(l) Selection of voltage program for source 11.2.3.linear variationisoadiabatic changesigmoid variationV(l) V(l) (V(l) V(l))*(m 1 sin 2nm)1/2moMoM ZnMcubic spline representation of arbitrary function froma table read from TAPE104.1KURVE(2) Selection of voltage program for source 2(as above for KURVE(l))NTlNumber of table entries forNT2Number of table entries forv l)v 2)program (ifNonepresent)programVlTABL(NTl,5) The V(l) program consists of NTl entires.The five elements of each entry are the nurn number mT,the relative amplitude at that turn (O VlTABL(mT,2) 0,and the three cubic spline coefficients for funcEionv l) in the neighborhood of mT:v l)m V(l) (V(l)o.MV l»1x [VlTABL(mT,2)]None VlTABL(mT,3) x (m:-mT) VlTABL(inT,4) x (m-In.r) 2 VlTABL(mT,5) x (m-mT) 3(2)V2TABL(NT2,5) Like VlTABL above but for VNonem4.5Tracking Cycle ParametersParameters principally related t.o details of .the particle tracking process areread according to the NAMELIST !CYCLE! in SUBROUTINE CYCPROG which is invoked bythe data command C. The possible input parameters follow:NAMEFUNCTION OR DEFINITIONUNIT.DEFAULTMTURNSNumber of turns to be tracked with fixed parameters100MSTEPNumber of tracking steps (minimum) per synchrotronperiod. (Active only for SWITCH(8) .TRUE., acceleratedcalculation feature - use cautiously)100ACCEL Number of beam turns per tracking step (maximum).(Active for SWITCH(8) .TRUE., accelerated calcu ation).1.
10I TRAPIndicate a condition for which tracking should beTM-12740interrupted before MTURNS0:1:2:3:4:5:LGRTHMSelect difference equation used in tracking1:2:3:4:SWITCH(i)no traptrap on minimum bunch widthtrap on minimum bunch heighttrap for n ETRAP (tolerance n/n .01)trap for 8 PHISTRP (tolerance Ps . 005)trap for n 0 - transition crossingcompleteSi Ss'completecomplete3kinematics, a(local) to ( p/p) 2ap I ap( p/p)kinematics, ap I ap( p/p)kinematics, ap to order p/pActivates optional program feature for SWITCH(i) .TRUE.(1) Voltage program called each step(2) Phase change program called each step(3) Adjust RF voltage for constant bucket height(4) Adjust RF voltage for constant bucket area(5) Time cycle in second (Tl and T2); number of turnsMTURNS must be M(T2)(6) Include a second RF source(7) Collect a history of variable parameters andproperties(8) Accelerate calculation by taking several beam turns(9) Calculate mean and rms width of theta and Eprojections of distribution(10) Phase feedback from bunch centroid (infinite bandwidth)(11) Turn-by-turn FFT of azimuthal distribution(12) Run RF source 2 at frequency determined by momentumdifference p/p DELTA2 from p(13) Use completely linear RF wavefo LSE.FALSE.FALSE.SDECRIf SWITCH(4) .TRUE. multiply bucket area by SDECR eachturn1.HDECRIf SWITCH(3) .TRUE. multiply bucket height by HDECR eachturn1.ETATRPIf ITRAP 3, stop tracking when ETA ETATRP.001PHISTRPI f ITRAP 4, stop tracking whenPHI SLIMIf SWITCH(3) or SWITCH(4) .TRUE. do not allow voltagereduction such that lsin P j PHISLIMsIsin Ps I PHISTRP.95.95
TM-12741 14.6Parameters Defining the Initial DistributionThe distribution of the particles to be tracked is established according toparameters read in according to NAMELIST IPOPLBI from the USBROUTINE POPUL8 which isinvoked by the data command P. The available parameters follow:NAMEFUNCTION OR DEFINITIONKINDChooses a type or mode distributionUNITDEFAULT11: Rectangular outline NTH by NE points, limited by THMIN,THMAX: TMIN, TMAX2: Uniform rectnagular grid NTH by NE, same limits3: Random uniform in THETA, Gaussian in E where rectangularlimits like KIND !; E limits are 204: Random uniform distribution of NTHxNE points withinrectangular limits like KIND l5: Gaussian in THETA, uniform random in E where rectangularlimits are extremes in E and 2a in THETA6 RecLanguLar rascer, regular In THETA, E spacing of rasterlines averages to Gaussian7: Matched bunch outline of NPOINT particles, SBNCH eVs8: Uniform, equispaced, matched bunch of appx. NTHxNEparticles9: Uniform,.random, matched of NPOINT particles10: Bi-Gaussian, NPOINT particles 95% within contour ofSBNCH eVs11: NPOINT uniformly spaced particles on flow lines justabove and below bucket boundary (not heavily used)12: Random uniform in E, parabolic in THETA13: Parabolic matched distribution of NPOINT particles,SBNCH eVsTHMINLower THETA limit on rectangular distributionradNoneTH MAXUpper THETA limit on rectangular distributionradNoneTMINLower kinetic energy limit on rectangular distributionMeVNoneTMAXUpper kinetic energy limit on rectangular distributionMeVNoneNTHNumber of grid points in THETA direction2NENumber of grid points in E direction2SBNCHArea of a matched or bunchIPOPSpecifies which RF source to use in matching1:2:3:THO FFeVs1.RF source 1 (eV 1 , w1 , h 1 )RF source 2 (eV 2 , w2 , h2 )Sum of sources 1 and 2Amount to displace distribution in THETA directionrad0.
12TM-1274NPOINTNumber of particles in matched distribution(except kind 8 where NTH and NE are used)PARTIONIf PARTION .TRUE., the distribution is partitioned irclasses which may be plotted with distinct symbols4.71.FALSE.Parameters Controlling Graphical OutputThe quantity and form of graphical output is determined by parameters read byNAMELIST !GRAPH! residing in SUBROUTINE GRAFSET.This subroutine is called by the0 data coIIlllland.NAMEFUNCTION OR DEFINITIONNPLOTNumber of turns between plots of distributionNoneIOPTSelects contour to appear on distribution plotNone0:1:2:3:4:5:IREFno contourbucket created bybucket created bybucket created bya matched contourSBNCHa matched contourUNITRF source 1RF source 2sum of source 1 and 2having the area of the initial bunch,containing 95% of the distributionSelects the RF source used in the matching for IOPT 4or 51:2:3:DEFAULT1RF source 1RF source 2sum of RF sources 1 and 2THPMINLower THETA limit for plotsradNoneTHPMAXUpper THETA limt for plotsradNoneDEPMINLower energy limit for plots in MeV from EsMeVNoneDEPMAXMeVNoneTITLEUpper energy limit for plots in MeV from EsIf TITLE .TRUE., read a title for graphical outputfollowing the NAMILIST [GRAPH! dataPLTSW(i)Plot options, PLTSW(i) .TRUE. activates feature.FALSE.FALSE.(1) Plot each class in the distribution (see PARTION indistribution parameters) with its own symbol andinterconned points of the same class(2) Like (1) above except omit interconnecting lines(3) Interconnect all points (regardless of PARTION)(4) Suppress captions, axis labels, etc.(5) Suppress symbols at individual points for PSTSW(i) .TRUE.(6) Replace scatter plot of distribution by six linesindicating mean and sigma for THETA and E projections(7) Suppre&s THETA histogram(8) Suppress E histogramDTHCURVMove contour on plot in THETA direction0.DECURVMove contour on plot in E direction0.
134.8TM-1274Selection of Parameters to be Recorded on History PlotsWhen SWITCH(7) in the cycle data is set, CYCPROG records a selection ofinteresting parameters for each tracking step on TAPE9.Any parameter can be plottedagainst any other by giving a pair of indices for each plot.The indices are elementsof the array NPLT(2,28) read in by NAMELIST IHISTRYj from subroutine HISTORY.subroutine is invoked by the H data command.ThisUp to 28 index pairs may be given.Theparameters are labelled by the index values as follows:NAMEFUNCTION OR DEFINITIONUNITNPLT(i,j) i 1, 2; j 1, 281:2:turn numbereVi, peak acceleration from source 1 per turn3:eV 24: 5:6: :25:26:27:28:1,,,peak acceleration from source 2 per turnphase of RF source 1phase of RF source 2 relative to source 12p', change of synchronous momentum8AV mean THETA of distributionEAV mean total energy of distribution09, rms width of distributionOE rms height of distributionEs, synchronous total energy Ps synchronous phase for source 1EL longitudinal emittance of distribution (onlycorrect for uncorrelated distributions as it nowstands; should be changed to be more useful)es, number of synchrotron periods per totation periodSB, bucket areaRB, bucket heightn, l/y - 1/y .ACCEL, the number of beam turns per tracking step inan accelerated calculation 8, the THETA of a particle which started from 0, Es;a measure of numerical error E, the E coordinate of the particle which started atthe synchronous point O, EsT, the beam circulation periodt, the integrated timeFamp 1 ; Fourier amplitude selected by NClFamp , Fourier amplitude selected by NC22Faze 1 , the phase of FampFaze2 , the phase of FampSpare 1 , open to plot. anything of interestSpare2 , adeVseVsMeVradMeVssradrad
14.94TM-1274Parameters Controlling FFTThe routine FOURFIT is invoked by the data command F.This routine will readdata according to NAMELIST IFFTI and analyze the current distribution using these data.If SWITCH(ll) is set in the cycle data, the same analysis will be carried out everytracking step and two selected amplitudes will be stored on the history tape.Theparameters are:NAMEFUNCTION OF DEFINITIONNBINFFTNumber of bins in representation of THETA distribution256NClNumber of amplitude and phase selected for historyh1NC2Number of amplitude and phase selected for history plotNone5.UNITDEFAULTPROGRAM AVAILABILITYThe program is available from the author as a CDC UPDATE source file in Fortran V.Because continual additions are made to the code there will be no guarantee that thematerial provided will be entirely consistent with this writeup, but basic featureshave been fairly stable for over a year.There are several hundred lines of comment inthe code which would help some one wishing to change it for special purposes.other documentation is available except for various examples of application.NoSee forinstance IEEE Trans. on Nucl. Sci .::v NS-30 114 p. 2627-2629 ("RF Exercises Associatedwith Acceleration of Intense Antiproton Bunches at Fermilab", J.E. Griffin, J.A. MacLachlan,Z.B. Qian), thesam .J ·2630-2633 ("Time and Momentum Exchange for Production andCollection of Intense Antiproton Beams at Fermilab", J.E. Griffin, J.A. 'MacLachlan,A.G. Ruggiero, K. Takayama), or Design Report Tevatron I Project, Sections 4.9-4.10and Ch 6.
5 TM-1274 If one synchronizes the RF frequency on each turn so that w hrl - hrl n o,n s,n then Pi,n l cp, 2Tih rls,n 1. ,n Q i,n l To simplify calculations in which there is more than one RF frequency present it is convenient to replace the RF phase variable by an azimuthal angle 8 measured from the RF source location in a direction opposite to the direction of particle circulation.
