TTRI-BN-15-03: Linac Envelope Optics
Rick Baartman, TRIUMFMarch, 2015
Abstract
I develop the formalism that allows calculation of beam envelopesthrough a linear accelerator given its on-axis electric field. Space chargecan naturally be added using Sacherer[1] formalism. A complicatingfeature is that the reference particle’s energy-time coordinates are notknown a priori. Since first order matrix formalism applies to deviationsfrom the reference particle, this means the reference particle’s time andenergy must be calculated simultaneously with the beam envelope andtransfer matrix. The code
TRANSOPTR [2] is used to track envelopes forgeneral elements whose infinitesimal transfer matrices are known, and inthe presence of space charge. Incorporation of the linac algorithm into
TRANSOPTR is described, and some examples given.
For beamline design, matching, etc., even in the case of strong space chargeforces, it is often sufficient to calculate the evolution of size parameters such asthe rms width and bunch length. The full set of such parameters is the6-dimensional “ σ -matrix”; x, P x , y, P y , z, P z . Conventionally, the momenta aredivided by the reference particle’s momentum, making them x (cid:48) , y (cid:48) , ∆ P/P .Parameter 5 is technically time, but converted to distance (bunch length) by Throughout this note, primes denote d/ds a r X i v : . [ phy s i c s . acc - ph ] S e p ultiplying by the speed β c of the reference particle. Similarly, parameter 6 is ∆ E , converted to P z by dividing by β c . The σ -matrix has 21 parameters, as it issymmetric: 6 rms sizes on the diagonal, 15 correlation parameters.The mathematical fomalism for this technique, including space charge, wasestablished by Frank Sacherer[1]. The space charge forces depend crucially uponthe bunch dimensions in configuration space, so it is important that these betracked. In other words, it is insufficient to use a formalism that first integratesthe equations of motion to derive the transfer matrices, and apply space chargeeffects afterwards. Some implementations such as TRANSPORT and
TRACE3D divide standard elements into (hopefully sufficiently short) segments interleavedwith space charge “lenses”. This is crude, approximate and non-adaptive.
TRANSOPTR [2] uses the envelope formalism, but did not until now include thecase of beam high intensity bunches being accelerated with a RF accelerator.This case is of particular importance for modeling the elecron linear accelerator.Short accelerator gaps can often be sufficiently well modelled as infinitesimal i.e.as a thin longitudinal lens with also a transverse thin lens focusing component.But this is insufficient for extended RF devices and as well misses time-of-flighteffects when the change in velocity of the reference particle is significant.Extended DC longitudinal fields have already been added to
TRANSOPTR [3]; thepresent study aims at adding the AC case. In the DC case, the reference particle’slongitudinal coordinates are inferred from the local potential. Not so the ACcase; the time and energy versus longitudinal location must be tracked in aseparate calculation and then the equations of first order deviations applied tofind the transfer matrices and local envelope optics.
This follows closely from the previous note on DC axially-symmetriclongitudinal fields[3]. 2 .1 Hamiltonian
With the distance along the reference trajectory s as the independent variable, theHamiltonian is H ( x, P x , y, P y , t, E ; s ) = − qA s − (cid:115)(cid:18) E − q Φ c (cid:19) − m c − ( P x − qA x ) − ( P y − qA y ) (1)The case of RF axially-symmetric electric field can be handled entirely with noelectric potential ( Φ = 0 ), and time-varying vector potential. This has beenpresented a number of times in the past (e.g. Chambers[4]), but we are interestedin the following more experimentally-useful case: The electric field along theaxis E ( s ) has been measured and is therefore known, and the geometry is exactlyaxially symmetric. Rob Ryne [5] has treated this case, and we use his vectorpotential (cid:126)A ( x, y, s, t ) directly. A x = E (cid:48) ( s )2 sin( ωt + θ ) ω x (2) A y = E (cid:48) ( s )2 sin( ωt + θ ) ω y (3) A s = (cid:32) −E ( s ) + x + y (cid:34) E (cid:48)(cid:48) ( s ) + ω c E ( s ) (cid:35)(cid:33) sin( ωt + θ ) ω (4)This is Coulomb/Lorenz gauge, satisfies Maxwell equations to second order intransverse coordinates, gives correct on-axis (cid:126) E = − ∂ (cid:126)A/∂t = E cos( ωt + θ ) .A priori, we do not know the reference particle’s energy and time coordinates.We need these in order to expand about them. They can be found from theequations of motion for x = y = P x = P y = 0 : dEds = ∂H∂t = q E cos( ωt + θ ) (5) dtds = − ∂H∂E = Ec P = 1 β c (6)(From here on, I drop the subscript: β and γ are implicitly assumed to be therelativistic parameters of the reference particle.)These are solved first and give the functions E ( s ) and t ( s ) about which t and E are expanded: E = E + ∆ E , t = t + ∆ t . So we transform the canonical3ariables t and − E to (∆ t, − ∆ E ) , using as generating function G = − (cid:32) t − (cid:90) dsβ ( s ) c (cid:33) (∆ E + E ) (7)(Check: ∂G∂t = − E , ∂G∂ ( − ∆ E ) = ∆ t .) The Hamiltonian gets the added terms ∂G∂s = ∆ E + E ( s ) β ( s ) c − ∆ tE (cid:48) ( s ) . Then expanding the square root, we get: H ∆ t = (cid:32) E βc − P (cid:33) − qA s − ∆ tE (cid:48) ( s )+ (∆ E ) β γ mc + ( P x − qA x ) + ( P y − qA y ) P (8)In expanding P x − qA x , P y − qA y , the time dependence disappears because it ishigher order: ( P x − qA x ) = P x − q E (cid:48) sin( ωt + θ ) ω xP x + (cid:32) q E (cid:48) ωt + θ ) ω (cid:33) x , (9)and similary for y . The term linear in ∆ t in the expansion of A s about t cancelsthe − ∆ tE (cid:48) ( s ) term, as it should but there is a remaining term quadratic in ∆ t ,the bunching effect. This leaves − qA s − ∆ tE (cid:48) ( s ) = q E sin( ωt + θ ) ω (cid:32) − ω (∆ t ) (cid:33) − r q (cid:32) E (cid:48)(cid:48) + ω c E (cid:33) sin( ωt + θ ) ω (10)Notice the first term here and the first term in eqn. 8 depend only on theindependent variable and not on the 6 dependent ones. Thus these do not affectthe equations of motion and we ignore them. We have: H ∆ t = − q E ω T (∆ t ) + (∆ E ) β γ mc − r q (cid:32) E (cid:48)(cid:48) + ω c E (cid:33) T + P x P − q E (cid:48) T xP x P + (cid:32) q E (cid:48) T (cid:33) x P + P y P − q E (cid:48) T yP y P + (cid:32) q E (cid:48) T (cid:33) y P (11)4e defined here T ( s ) = sin[ ωt ( s ) + θ ] /ω to clean up the notation a bit.Finally, we wish to transform from (∆ t, − ∆ E ) to ( z, P z ) = ( − βc ∆ t, ∆ E/ ( βc )) .(The reason for the sign change is as follows: an early arrival implies ∆ t < , butthis means the particle is ahead so z > .) The generating function is G = − βc ∆ tP z (12)(Check: ∂G∂ ∆ t = − ∆ E , ∂G∂ ( P z ) = z .) The term to be added to the Hamiltonian is ∂G∂s = β (cid:48) β zP z = γ (cid:48) β γ zP z = q E CβcP γ zP z , where C ≡ cos( ωt + θ ) . H z = P x P − q E (cid:48) T xP x P + P (cid:32) q E (cid:48) T (cid:33) − T (cid:32) q E (cid:48)(cid:48) + ω c q E (cid:33) x P y P − q E (cid:48) T yP y P + P (cid:32) q E (cid:48) T (cid:33) − T (cid:32) q E (cid:48)(cid:48) + ω c q E (cid:33) y P z γ P + 2 q E Cβc zP z γ P − q E β c ω T z (13) Ryne[5] has a transformation that gets rid of the second derivative of the on-axiselecric field. It’s complicated. At the same time he transforms away the adiabaticdamping; it’s a neat and didactic trick but not strictly necessary forcomputational purposes. It is simple to just use P x,y,z directly and then justrescale by final P at the end.But there’s an easy way to get rid of the second derivative: it turns out that thevector potential can be simplified if we use a different Gauge.I propose the following function Ψ( x, y, s, t ) = − E (cid:48) ωt + θ ) ω x + y (14)5dd the gradient of this function to the previous vector potential (2,3,4). Thiszeroes both A x and A y , leaving A s = −E ( s ) (cid:32) − ω c x + y (cid:33) sin( ωt + θ ) ω (15)This is considerably simpler, but now there is a scalar potential: Φ = − ∂ Ψ ∂t = E (cid:48) cos( ωt + θ ) x + y (16)Now if we expand the Hamiltonian, we get a different result: H z = P x P + P y P + q βc (cid:32) E (cid:48) C − E S ωβc (cid:33) r P z γ P + 2 q E Cβc zP z γ P − q E ωSβ c z (17)( C ≡ cos( ωt ( s ) + θ ) , S ≡ sin( ωt ( s ) + θ ) )This is not only much simpler than eqn. 13 ( P x and P y have their usualdefinitions, no transverse cross terms, no E (cid:48)(cid:48) ), but has nice intuitive explanationsfor the individual terms. (1) The factor in parentheses is precisely the integrandof eqn. 4 of a short note I wrote in 1985[6] explaining the derivation of the focalpower of an RF gap, e.g. a buncher. (2) Taking the limit as ω → reproducesprecisely the Hamiltonian of the DC accelerator I derived in 2010[3]. Note thatin that case, E (cid:48) = − φ (cid:48)(cid:48) . A convenient and useful way of representing the equations of motion through theoptical element is the so-called infinitesimal transfer matrix approach[1]. Theinfinitesimal transfer matrix F ( s ) is defined as ( T − I ) /ds where T is thetransfer matrix from s to s + ds and I is the identity matrix. With this definition6ne has for individual particles d x ds = F x , where x ≡ xP x yP y zP z . (18)Beams of particles are conveniently represented by the so-called σ -matrix; theelements of which represent second order moments of the beam[1]. The σ -matrixand the transfer matrix M transform through the system according to theequations dσds = F σ + σF T , (19) dMds = F M. (20)where F T is the transpose of F .Now that the Hamiltonian for linear motion (eqn. 17) has been obtained, it is asimple matter to find the infinitesimal transfer matrix F . Writing the equations ofmotion ( x (cid:48) = ∂H/∂P x , P (cid:48) x = − ∂H/∂x , etc.) in the form of Eqn. 18, thefollowing F -matrix is found for the axially symmetric linear accelerator: F = P A ( s ) 0 0 0 0 00 0 0 P A ( s ) 0 0 00 0 0 0 β (cid:48) β γ P B ( s ) − β (cid:48) β . (21)where we have defined: A ( s ) = − q βc (cid:32) E (cid:48) C − E S ωβc (cid:33) , (22) B ( s ) = q E ωSβ c . (23)7 .4 Space Charge Space charge forces depend recursively upon the σ -matrix elements, and aresimply added to the focusing elements F n, m − | n,m =1 , , of the element’sinfinitesimal transfer matrix such as eqn. 21 above. This technique is used in thecode TRANSOPTR , as described by de Jong[2]. The resulting equations can onlybe solved numerically.The given references[1, 2] treat space charge in the non-relativistic regime, so itwas not obvious that
TRANSOPTR was correct in the relativistic regime. Thereare two effects that need to be considered to generalize the equations: the spacecharge magnetic field, and bunch length contraction. For detailed derivations, theinterested reader is referred to the Ph.D. thesis of Fubiani[7]. The first effectrequires dividing the space charge force by γ . The second requires that theCarlson elliptic integrals’ arguments be modified. The 3 integrals for the 3 majoraxes are R D ( u, v, w ) = 32 (cid:90) ∞ dt ( t + w ) (cid:113) ( t + u )( t + v )( t + w ) where ( u, v, w ) are ( σ , σ , σ ) for the x -axis, ( σ , σ , σ ) for the y -axis, ( σ , σ , σ ) for the z -axis. To be relativistically correct, σ must be replacedby γ σ . (See Fubiani[7], Appendix J.)An interesting limit is the long bunch, since this can be approximated as acontinuous beam with current I . First of all, it is clear that for this limit to apply,bunch length (cid:29) transverse size is not a necessary condition. Rather, γ timesbunch length (cid:29) transverse size. This means that for example a 1 mm long by1 mm wide electron bunch is already well into the long-bunch regime withenergy of 10 MeV.Secondly, in the long bunch regime, the Carlson integrals governing transversespace charge are ∝ ( γ √ σ ) − , or the inverse bunch length augmented by thefactor γ . For this reason, the constant governing space charge force in TRANSOPTR in the unbunched case is divided by an extra factor of γ . Note that this assumes the bunch axes are aligned with the beam direction ( s ) and the chosentransverse orientation. If this is not the case, a rotation must be applied. Implementation into
TRANSOPTR In TRANSOPTR , the momentum is dimensionless and in absence of accelerationcorresponds to angles ( x (cid:48) , y (cid:48) , z (cid:48) ) . When acceleration occurs, angles are not scaledcanonical momenta. The simplest implementation is to scale ( P x , P y , P z ) byinitial total momentum P i . Then after integration is complete, the momenta canbe converted back into angles by multiplying by P i /P f . In this case, F (eqn. 21)becomes: F = P i P A ( s ) P i P i P A ( s ) P i β (cid:48) β P i γ P B ( s ) P i − β (cid:48) β . (24)In order to find the σ -matrix at any point s of a beamline, it is sufficient tocalculate the transfer matrix to that point. This technique is not useful whenspace charge forces are non-negligible since then the optics themselves dependupon the σ -matrix. In the past space charge has been added to legacy codes suchas TRANSPORT by inserting a sufficient quantity of defocussing lenses whosestrength depends upon the local beam size. This is inefficient at best: it is thecrudest sort of integration. The best way to do this is as implemented in
TRANSOPTR : the exact Sacherer 6D envelope equations (19) are integrated witha higher order integrator such as the Runge-Kutta.In many situations, one needs only the σ -matrix and so it is sufficient to solve the equations of eqns. 19. However, in general one would like to have the abilityto fit certain optical characteristics, such as point-to-point imaging. It is thereforeoften useful to also know the transfer matrix to any point. This is in fact the incoherent transfer matrix. So in order to cover all cases,
TRANSOPTR solves In fact, since this matrix is symmetric, one can get away with solving only equations. Onecan whittle this down further for uncoupled cases; in fact first order or second order equationsare often all that is needed for a straight beamline with no coupling and no dispersion. These areof course the Kapchinsky-Vladimirsky equations. equations. Internally, the code uses asingle by matrix.Without time-varying fields, the total momentum of the reference particle isalways known at any s . In a linac, this is no longer true. Thus, theimplementation needs to track not only the transfer and σ matrix elements, butalso time and energy. These can be thought of as the 5th and 6th “zero-order”coordinates, with respect to which the fifth and sixth phase space coordinates aremeasured. Analogously, in the most general cases, the first through fourth“zero-order” coordinates may also not be known. This could happen for exampleif given only the field map of a dipole magnet. Thus there could be 6“zero-order” coordinates as a function of s , and the usual 6 phase spacecoordinates are differentials with respect to these.I therefore have added a row to the 72-element matrix, making it dimensioned13-by-6, and there are 78 ODEs. This would allow at a later date addingelements whose reference trajectory is not known a priori. Here follows the f77 subroutine that takes the 13-by-6 matrix SX and calculatesthe 13-by-6 matrix DSX , which are the derivatives of SX with respect to Z = s .Note that SX(13,5) is not time t directly but is ct , and SX(13,6) is notenergy directly but is γ − . The F matrix is in fact F directly. SUBROUTINE SCLINAC(Z,SX,DSX)COMMON/PRINT/IPRINT,IQ1,JQ1,IQ2,JQ2,IQ3,JQ3,IQ4,JQ4COMMON/VELEM/VELM(21,9999),NELMDIMENSION SX(13,6),DSX(13,6)DIMENSION TEMP(6,6),ROT3(3,3),ROT6(6,6),SIG(3,3)DIMENSION WRK(3),EVAL(3),RF(6,6),FT(6,6)LOGICAL STRATEcommon/splind/klob,kloeCOMMON/FMATRX/F(6,6)COMMON/SCPARM /KSC,ISCCOMMON/CALLS/NSC OMMON/MOM/P,BRHO,pMASS,ENERGK,GSQ,ENERGKi,charge,bnchargeCOMMON/axezrf/omoc,phase,etot,$ ZTBLE(9999),EZTBL(9999),SPCFE(9999),zse,iez,iaxezrfCOMMON/axez2/pratioREAL KSCdo i=1,4dsx(13,i)=0.enddoNSC=NSC+1pratio=1.gamm1=sx(13,6) ! gamm1=gamma-1energk=gamm1*PMASSif (energk.le.0.)thenwrite(6,*)’*********BEAM DECELERATED TO A STOP’write(6,*)’Check parameter values or limits.’call erroutendifc call spline routine to find the on-axis electric field and derivativescall spling(kloe,ZTBLE,EZTBL,SPCFE,IEZ,z-zse,EZ,EZp,EZpp)ccc=cos(sx(13,5)*omoc+phase)sss=sin(sx(13,5)*omoc+phase)c omoc is omega/cdsx(13,6)=EZ*ccc/PMASS/100. ! denergy/ds ; EZ is in MeV/m, already contains CHARGEgamma=1.+gamm1ETA=SQRT(gamm1*(gamm1+2.))BRHO=ETA*PMASS/CHARGE*3.3356397BRHO=BRHO*.001P=ETA*PMASS !this P is actually Pcgsq=GAMMA**2BETA=eta/gammadsx(13,5)=1./beta ! dctime/ds=1./betagamm1old=etot/pmassETAold=SQRT(gamm1old*(gamm1old+2.))Pold=ETAold*PMASSpratio=etaold/etac re-calculate space charge parameter as energy changes.c (898755.2 is 1/(4pi epsilon_0) ˜ cˆ2/1.e7 in some units)KSC=BNCHARGE*CHARGE*898755.2/(ETA**2*PMASS)/pratio (1,2)=pratioF(3,4)=pratioF(2,1)=(EZ*sss*omoc-EZP*ccc/beta)/(200.*Pold)F(4,3)=F(2,1)c bpob = beta’/betabpob=dsx(13,6)/(gamma*eta**2)F(5,5)=bpobF(5,6)=pratio/gsqF(6,5)=EZ*sss*omoc/(100.*Pold*beta**2)F(6,6)=-bpobIF (KSC.NE.0.)THENc *******************************************************************c START SPACE CHARGE SECTION:c *******************************************************************STRATE=.TRUE.c SIG is real-space alone (3D)DO 5 I=1,2DO 5 J=1,25 SIG(I,J)=SX(2*I-1,2*J-1)gamma=sqrt(gsq)C QUICK AND DIRTY Lorentz transform.SIG(1,3)=GAMMA*SX(1,5)SIG(2,3)=GAMMA*SX(3,5)SIG(3,1)=SIG(1,3)SIG(3,2)=SIG(2,3)SIG(3,3)=GSQ *SX(5,5)c write(36,*)’gsq’,gsq,sig(3,3)c *******************************************************************c Test whether bunch axes are along coordinate axesc *******************************************************************IF (SIG(1,2).NE.0..OR., SIG(1,3).NE.0..OR., SIG(2,3).NE.0.)THENCALL PRESET(ROT6,36,0.)STRATE=.FALSE.c Find rotation requiredcall jacobi(sig,3,3,eval,rot3,nrot)c Fill 6D rotation matrix O 55 I=1,3DO 55 J=1,3ROT6(2*I-1,2*J-1)=ROT3(I,J)55 ROT6(2*I,2*J)=ROT3(I,J)ELSEc *******************************************************************c Bunch is along axesc *******************************************************************DO 52 I=1,352 EVAL(I)=SIG(I,I)ENDIFdo i=1,3if (eval(i).le.1.e-6)theneval(i)=1.e-6endifenddoc RD is the Carlson elliptic integral R_D.c Scale the arguments. Carlson routine is better if args are order 1.scalrd=eval(1)if(eval(2).gt.scalrd)scalrd=eval(2)if(eval(3).gt.scalrd)scalrd=eval(3)rscalrd=scalrd**(-1.5)EL1=rscalrd*RD(EVAL(2)/scalrd,EVAL(3)/scalrd,EVAL(1)/scalrd)EL2=rscalrd*RD(EVAL(3)/scalrd,EVAL(1)/scalrd,EVAL(2)/scalrd)EL3=rscalrd*RD(EVAL(1)/scalrd,EVAL(2)/scalrd,EVAL(3)/scalrd)DO 70 I=1,6DO 70 J=1,670 RF(I,J)=0.RF(2,1)=KSC*EL1RF(4,3)=KSC*EL2RF(6,5)=KSC*EL3*gsqc *******************************************************************c no rotation requiredc *******************************************************************60 IF (STRATE)THENDO 40 I=1,6DO 40 J=1,640 FT(I,J)=F(I,J)+RF(I,J) LSEc *******************************************************************c rotate to lab framec *******************************************************************DO 15 I=1,6DO 15 J=1,6TEMP(I,J)=0.DO 15 M=1,615 TEMP(I,J)=TEMP(I,J)+RF(I,M)*ROT6(J,M)DO 20 I=1,6DO 20 J=1,6FT(I,J)=F(I,J)DO 20 M=1,620 FT(I,J)=FT(I,J)+ROT6(I,M)*TEMP(M,J)ENDIFc *******************************************************************c END SPACE CHARGE SECTIONc *******************************************************************c START NO SPACE CHARGE SECTIONc *******************************************************************ELSE !ie. KSC=0.DO 23 I=1,6DO 23 J=1,623 FT(I,J)=F(I,J)ENDIFc *******************************************************************c END NO SPACE CHARGE SECTIONc *******************************************************************c MATRIX OPERATIONS TO OBTAIN DIFF. EQUATIONSc *******************************************************************DO 30 I=1,6DO 30 J=1,IDSX(I,J)=0.DO 25 M=1,625 DSX(I,J)=DSX(I,J)+FT(I,M)*SX(M,J)+SX(I,M)*FT(J,M)30 DSX(J,I)=DSX(I,J)DO 26 I=7,12DO 26 J=1,6 SX(I,J)=0.DO 26 M=1,626 DSX(I,J)=DSX(I,J)+FT(I-6,M)*SX(M+6,J)RETURNEND
The TRIUMF injector electron linac, EINJ, takes bunches from energy of300 keV to ∼ MeV if properly phased and the peak gradient is MV/m.Below is example for phase θ = 0 at the start of the calculation. S t R t e n t r (c) R. Baartman 05/26/15 E I N J e x it Red is the 2rms transverse size, and green is the 2rms longitudinal (bunchlength). The input bunch parameters are somewhat arbitrary, roughly thecondition for a minimum beam size at exit. This particular case has zero bunchcharge. 15n this second example,
TRANSOPTR is instructed to fit the matrix element tozero. This makes energy insensitive to input phase, thus finding the peak energygain phase. This phase turns out to be θ = − . ◦ . S t R t e n t r (c) R. Baartman 05/26/15 E I N J e x it In the third example, bunch charge has been raised to pC.16 S t R t e n t r (c) R. Baartman 05/26/15 E I N J e x it Each calculation above takes roughly 400 Runge-Kutta steps for 2400 calls to the
SCLINAC routine. This gives 5-figure accuracy to the transfer matrix and the σ -matrix, and is easily enough for describing reality considering that the on-axisfield is only known to 2 or 3 significant figures. The extra accuracy is useful,however for fitting matrix or beam matching, which is done with a downhillsimplex method, or simulated annealing for cases of more than 3 fittingparameters.On my unremarkable, circa 2006 Intel desktop, each run through the linac takesabout 17 milliseconds with zero bunch charge and 25 milliseconds with spacecharge. The difference is due to the Carlson elliptic integrals needed for thespace charge case.On a typical optics matching case, one varies 2 solenoids, the buncher amplitude,17nd the linac phase, to minimize the bunch size and energy spread at the linacoutput. A calculation with such a fit requires typically a half million total calls to SC (the space charge routine for no-linac case) and SCLINAC , and so takes about5 seconds CPU time. The result is shown below. The bunch charge is pC.Rather extravagantly, each calculation starts from the cathode whereas it wouldhave been more efficient to store the beam parameter set at the buncher entranceand start it from there. The DC acceleration to 300 keV from the cathode isdescribed in reference [3].The Buncher itself, located at s = 85 cm, is calculated as just another linac,phased to give no energy gain. 18 S t R t e n t r (c) R. Baartman 05/26/15 A cc l e x it S O L g e n t r B N C H e x it S O L S O L M a t c e n t r E I N J e x it M db0 e-, from restx,y_2rms/mmz_2rms/mmdp/p_2rms/%Energy/MeV The code describing the whole transport from cathode to linac exit, is shownbelow.
SUBROUTINE TSYSTEMCOMMON /BLOC1/c0,Ebun,phbun,c1,c2,phase,grad,birkCOMMON/MOM/P,BRHO,pMASS,ENERGK,GSQ,ENERGKi,charge,currentCOMMON/SS/SX(13,6)nss=2 sollen=6.32 this was with old clampssollen=5.64sollenrot=7.29solt=c0*140.0/10000.sol1=c1*116.9/10000.sol2=c2*116.9/10000.c following are from autocad dwg.atgunsol=42.193atdb=56.9atrfgap=85.642atdb2=105.7atsol1=122.777atsol2=176.577atmatch=237.sap=2.5/2.*2.54sss=sollen/nssc since actually a solenoid’s aberration is unrelated to its length,wt=ww/sqrt(float(nss))acclen=12.0call axez(59,121,-100.,acclen,10,1,2)call drs(atgunsol-acclen-sollen/2.)rho=100.*brho/soltth= sss/rho/2.*57.29578do i=1,nssif(i.eq.(nss+1)/2)thencall solenoidn(solt,sss,sap,wt,’SOLg’)elsecall solenoidn(solt,sss,sap,wt,’ .’)endifcall rotate(-th,’ .’)enddoth= sollenrot/rho/2.*57.29578call rotate(th,’ .’)if(iprint.ne.0)write(6,*)’Sol.B/1gauss,rot.= ’,solt*10000.,th*nsscall drs(atdb-atgunsol-sollen/2.)call drs(atrfgap-atdb-10.) !-ccycall linacn(62,21,Ebun,20.,1.3e9,phbun,20,’BNCH’) !BuncherEHsollen=5.85 ollenrot=7.37sss=sollen/nsscall drs(atdb2-atrfgap-10.)call drs(atsol1-atdb2-sollen/2.)rho=100.*brho/sol1th= sss/rho/2.*57.29578do i=1,nssif(i.eq.(nss+1)/2)thencall solenoidn(sol1,sss,sap,wt,’SOL1’)elsecall solenoidn(sol1,sss,sap,wt,’ .’)endifcall rotate(-th,’ .’)enddoth= sollenrot/rho/2.*57.29578call rotate(th,’ .’)if(iprint.ne.0)write(6,*)’Sol.B/1gauss,rot.= ’,sol1*10000.,th*nsscall drs(atsol2-atsol1-sollen)rho=100.*brho/sol2th= sss/rho/2.*57.29578do i=1,nssif(i.eq.(nss+1)/2)thencall solenoidn(sol2,sss,sap,wt,’SOL2’)elsecall solenoidn(sol2,sss,sap,wt,’ .’)endifcall rotate(-th,’ .’)enddoth= sollenrot/rho/2.*57.29578call rotate(th,’ .’)if(iprint.ne.0)write(6,*)’Sol.B/1gauss,rot.= ’,sol2*10000.,th*nssd5=74.96-sollen/2.call drs(atmatch-atsol2-sollen/2.)call drift(0.,’Matc’)call drs(15.-1.023)call fit(1,1,1,0.,1.,1)call fit(1,5,5,0.4,8.,1)call drift(1.023,’ .’) x(13,5)=0.call linacn(71,447,grad,127.953,1.3e9,phase,100,’EINJ’)call drift(1.023,’ .’)call fit(1,1,1,0.,4.,1)call fit(1,2,2,0.,1000.,1)call fit(1,5,5,0.,70.,1)call fit(1,6,6,0.,1000.,1)enrg=sx(13,6)*0.511c we don’t want output energy to be too far off.call fitarb(10.02,enrg,100.,1)call drs(61.37)call drift(0.,’Mdb0’)call vective(1)returnend References [1] Frank J Sacherer. Rms envelope equations with space charge. Technicalreport, CERN DL/70-12, 1970.[2] MS De Jong and EA Heighway. A first order space charge option fortransoptr.
Nuclear Science, IEEE Transactions on , 30(4):2666–2668, 1983.[3] R. Baartman. Bunch dynamics through accelerator column. Technicalreport, Technical Report TRI-DN-10-01, TRIUMF, 2010.[4] EE Chambers. Particle motion in a standing wave linear accelerator.Technical report, DTIC Document, 1968.[5] RD Ryne. The linear map for an rf gap including acceleration. Technicalreport, LANL Report 836 R5 ST 2629, 1991.[6] R. Baartman. Buncher Beam Dynamics. Technical Report TRI-DN-85-37,TRIUMF, 1985.[7] Gwenael J Fubiani. Controlled electron injection into plasma acceleratorsand space charge estimates.