A program package for designing accelerating structures

U.S.S.R. Comput.Maths.Math.Phys., Printed in Great Britain

0041-5553/86 $lO.OO+O.OO 01988 Pergamon Journals Ltd.

~01.26,No.4,pp.159-164,1986

A PROGRAM PACKAGE FOR DESIGNING ACCELERATING STRUCTURES'

A.G. DAIKOVSKII,

YU.1. PORTUGALOV

and A.D. RYABOV

Naves of any type are calculated in periodic structures. In accordance with Floquet's representation, a problem is formulated for an element of periodicity and, when there is symmetry, for half of it. E-node tetragonal finite elements are used for the discrete approximation of The algebraic eigenvalue problem is solved the differential problem. using the subspace iteration method. The determination of eigen electromagnetic waves in a volume bounded by an ideally conducting metal surface consists of obtaining the eigenvalues and eigenfunctions of the homogeneous Maxwell equations. This problem has numerous applications, of which a fundamental one is the design of accelerating structures in charged-particle accelerators. There are a number of programs for solving the problem numerically /l-3/, and we can divide them into two categories. The first covers programs which use the net method (finite differences of finite elements). The second category is based on the particle domain method. The latter was used, for example, in /4-7/, where references to earlier publications can be found. 1. Formulation For fields which

of

the

problem.

change with time as e'"*, Naxwell's VxE=-ikH,

The boundary

condition

V.E=O.

ikE,

VXH-

on the metal

equations

have the form

V.H=O.

(1)

is EXv=O.

(2)

Here E,H are the complex amplitudes of the field, of periodic guiding structures we can seek the eigenwaves

k=olc,v is the normal. In the case in the Floquet representation

H(z, y* z)=e-‘@=“#(z,

E(s, y, z)-~-‘~*‘=EP(~, y, z),

y, z),

where 8, % are periodic (with the period L) complex vector-functions, and L is the geometric period of the structure along the s axis. From these relations there follow the equations which we call the conditions of quasiperiodicity:

E(z,y,z+L)=e-“E(~,y,z),

H(z,y,z+L)=e-‘“H(r,y,z).

(3)

We shall confine ourselves to real values of the parameter B,O
(4)

O~\V'+kpW-0 with boundary

conditions

on the metal

of the first kind wxv==o,

P.W-0

(5)

or of the second kind [VXW]Xv==O, where W denotes E or 11 repsectively in (5) and - again with the conditions of quasiperiodicty

W.v-0, (6), and in the case of a periodic

W(z, y, z+L)=e-“W(z,

y. 2).

(6) structure

(5)

The system of eigenfunctions of problems (4)-(7) is, as is well-known /9/, the direct sum of the vortex and potential subsystems. Therefore, at least for non-multiple eigenvalues,

or

*Zh.vychisl.l4at.mat.Piz.,26,8,1206-1214,1986

159

160 If one and the same eigenvalue corresponds to both the vortex and potential eigenfunctions, we are concerned with the invariant subspace drawn on these functions. In this case we can recommend proceeding from E-components to H-components or vice versa. The vortex function will at the same time correspond to one and the same electromagnetic field (E. H), and the potential subsystem spectrum changes. If the resonator has a plane of symmetry, it is sufficient to consider an elementary fragment of the domain, extending the solution synnnetrically or antisymmetrically through these planes to the whole domain. Boundary conditions of the first kind correspond to an antisymmetric extension, and those of the second kind correspond to a symmetric extension. We shall distinguish the cases for the periodic structure when we can choose the plane of symmetry :=consl for a fragment of periodicity and when this plane does not exist. Turning to (4)-(71, we see that W is determined to within a phase factor. This arbitrariness can be eliminated in the formulation of the problem for a symmetric period. If L/2
Im+W,=O, ”

ImW,=O,

Re&W,=O, Reg-WY-o, in the plane

z=O

Re

w,=o

and Im (eie12Wr) = 0, Re j&+

in the plane

(84

w=) = 0,

im(eief*W,)=

0,

1111 i&2.$

a Re(e@fl - W - 0 az It-’)

Wz) = 0,

Re (etR’dWz) = 0

z=L/2.

2. Solution of the algebraic eigenvalue problem. The numerical solution of problem (4)-(7) will involve its discretisation finite element method and its reduction to an algebraic system of the form

using the

Aw=l.Bw (10) A and B. The complexity of the latter with generally Hermitian, positive-definite matrices In problem (10) it makes sense not to is connected with the condition of quasiperiodicity. seek all the eigenvalues and eigenvectors, but only those which relate to the low-frequency The inverse iteration method is used below with the simultaneous part of the spectrum. iteration of several eigenvectors to solve this partial problem (see /lo-13/). Take For the real (nXn) matrices A and B the subspace iteration method is as follows. Then the following matrix X, m linearly independent vectors, and suppose they form the nXm holds in the k-th step of the iteration process: 1) the system AY,-&BY,==BX,_, is solved; 2) the matrices C,==Y,‘AY,, D,-Y,‘BY, are formed; 3) the minor eigenvalue problem C,Qp=D,Q&,. is solved; is calculated. 4) XF-Y,Q, approaches the invariant subspace drawn to m eigenvectors The subspace which, as k-m, with eigenvalues which are close to ho is separated at stage 1. Stages 2)-4) are the Rayleigh-Ritz procedure for obtaining the best approximations to these eigenvalues and eigenStage 1) is implemented in the form of an LU-expansion vectors from the iterated subspace. A--hOB. In Stage 3) we can use standard programs based on a QR-algorithm or on of the matrix Note that, by obtaining eigenvectors of the minor problem which are Jacobi's method /12/. D*-norm, we automatically obtain orthonormalised Rayleigh-Ritz orthonormalised in the and consequently, in the limit as k-t=, also the eigenvectors approximations, X,TBX,-Q,TY~‘BY*Q~-Q~TD,Q,=I. In the case of the Hermitian

matrices A,‘=A,,

We shall represent

A and A,‘=-A,,

B,A=A,+iA,,B=B,+iB,, B,‘=B,,

them in the form of synrnetric extended

we have

B,‘=-B,. real matrices:

It is obvious that for each eigenvalue ?, together with the eigenvector W-I& ~1' the vector is orthogonal to 1~'in the B-norm (w’Bw=O), is an eigenvector, i.e. each u=:l-u, uj', which eigenvalue is at least doubly degenerate. The properties of the matrices A and B enable us to implement stages l)-4) economically. Matrix X is a (LnXm, n)-matrix, Suppose X is again a matrix of m linear independent columns. We shall form the matrix IlXiXjl,which we shall use for the initial i.e. a (u. v)-matrix. Then (the k index is omitted) the following holds: approximation. is not solved, since 1') the system AY--haBY-BX is solved, and the system AY-?.JJF-BX y=Y (an economy) ; 2') the matrices C=lYiFj/‘A[jYiFll, D==~~YiP~~‘B~~YiP~~,are formed, at the same time

161

where C,=l"AI'.D,=Y'BI‘.C:=YAP, D,=Y'BY, the matrices C,, D,are symmetric symmetric (a saving in convolutions when forming matrices C and D); 3') the following minor problem is solved;

andC2,D,

are anti-

CQ=DQQ; it

is

easy

to see that the matrices C and D inherit the properties of matrices A and B, since the minor problem has the form

am) is the matrix of the eigenvalues, and Q=IIRiRI/is the matrix of the where Q,=diag (OS,..., eigenvectors of the minor problem; 4') X=(/YiY/jRis calculated (the saving here is that x is not calculated). The modification discussed in l'f-4'1 enables us to halve the volume of calculations and to economize an computer memory. 3. Axisymmetric structures. Suppose there is an axisymmetric structure and wdenotes e or h, where e=(e,,e,,.e,} and h=(h:. h*, &} are the field amplitudes obtained after separating the harmonic dependence with respect to cp,i.e. they are only dependent on the coordinates z and P. Then the wave Eq.(4) takes the form (V=(a/&,a/@))

(,llb) (lid

We will assume n#O since the case with n-0 reduces to a simpler scalar problem. We need system (11) in fairly general notation, more precisely - in projections on to the (z,, rp,rp) axes, rotated in the (2,~) plane by an angle Q with respect to the (a,Ptm) axes, such that w has the components (EL,a%,u;). It is easy to obtain the following system:

- QpVwl + ap 9 1 ----_*+--zc,= h ata P

(ll’a)

ap 2n p

6%

0,

-aP 2n m&=0 32s P ’

ll’c)

Note that ap/ax,=sina, dP/8x~=cesa. We shall write boundary conditions (5) and (6) for system (11'). Suppose fv,t) is a local mobile set of coordinates which is naturally connected with the boundary contour in the (3,P) plane. Then the conditions of (5) take the form

where I is the local radius of the contour curvature, whilst the I*+"sign corresponds to the convexity of the boundary to the side ofthenormal, and otherwise the 11-11 sign is used. Accordingly, the conditions (6) have the form

We recall that (12) and (13) relate to the boundary conditions on the metal, if, correspondingly, or w==h. In the planes of syannetryz=const either (12) or (13) are chosen as a w-e function of the synrnetryof the initial modes which in this case look particularly simple:

wml

26:4-K

162 on the axis

tw,~w~=o, n> f. The discretisation and reduction of problem (ll')-(16)to an algebraic problem is carried out in a standard way using the finite element method /la/. Suppose a section of the structure in the (z,p) plane is divided into finite elements with some numbering of the nodes and some set of basis functions {@(z,p)}.We shall represent w in the form i,i-l, Z,..., IV, of an interpolant:

We shall project the nodal values of 4 in all the internal grid nodes on to the (2,P? 9) axes, and at the nodes at the boundary of the metal - on to the direction (v,T,cp). Multiplying (11') bu Jliand integrating bj parts, at each internal node we will obtain an equation for w,j, wO*$ w,f:

17) where

Equations for (wol',u)tflwV‘ ) follow from the boundary node i from (11'): bi;sinacosa

2nbijsina

where

(19)

In (19) when integrating by the boundary of the element we can replace iYw,/& and aw#v by aw,/& and aw,/dvto within 0(&/r) and can use the conditions (12) and (13); 6c is the linear dimesnion of the element. In order finally to close the problem it remains to establish the connection between the projectionsw'at the internal nodes with the projections at the nodes at the boundary in Eqs.(l7) and the inverse connection in Eqs.(l8). Suppose f&i) is a fixed pair of indices which occur simultaneously in the sums (17) and (18), and suppose i is an internal node, and j is a boundary node. Bearing in mind that wi and WI are projected in different systems of coordinates in (17), (181, we make the following substitutions:

a

Fig.1

b Fig.2

163

Then the matrix elements of the connection of the components of tv'with IV'in (17) and \\'with using direct wP in (la! equal M,P and M,:P respectively. It is easy to obtain (M,,P)'=-M.,'Y' multiplication, which indicates the symmetry of the complete matrix.

Fig.3 Remark

1. Taking the conditions of quasiperiodicityof (7) for the arbitrary period, and

of (8) and (9) for the symmetric period instead of the conditions of symmetry of (14) and (15), we can similarly show the symmetry of the matrices of the discrete problem. 2. For the approximations in the neighbourhood of theconvex angle (the condition at the edge /4, 18/) to be correct, the latter is assumed to smooth out and locally to refine the net in order that we can disregard the term O(‘sc/r). Table 1

Table 2

bfv 6.7.10-' 1.2. io-’

3.1'10-* !.0.10-0 4.4.10-'

C’&I

1

1.3.10-5 2.5.10-a ;+:;g 7:*.10-:

t&x

1

4.0.10-3 i.%lO-' 1.7.w r,.ti. IO-' 2.1.10-'

bfr:

)

4.6.10-J 2.3.iW+ 5.2.10-' l&IO-' 6.9.Kw

&

G.?. lo-:

2.i;.Io-a

6.4.iO-' %.I.lO-' 8.5.10-t

Fig.1 shows the parametrisation of the acclerating structure half-cell with washers and diaphragms. The calculated lower dispersion branches are given in Fig.2 for the values of the parameters r,,rl,r,,r..r,,which equal 190, 137, 133, 35 and 17 millimetres respectively, and for t,,I,,t,,which equal 8,5,22 (a-30"). At the same time Fig.2a corresponds to waves with one variation in azimuth, and Fig.2b corresponds to waves with two variations. Fig.3 shows the isolines e,and & for 9=n!4 of the two lower branches of Fig.2a. The fields in the cells are reduced using quasiperiodicity with respect to the values in the first half-cell. Fig.4 shows the vector field (e.,e,)of a wave of the form 8=2z/3 of the EH, tx= for one of the versions of the structure.

4. ikief data on the program package. The above procedure is implemented in a program package with the general title PRUD-W. which consists of the following modules:

164 an octagonal isoparametric element net generator; a basic program system for obtaining eigenfrequenciesand fields; a module for calculating the secondary quantities necessary for applications - the ratio characteristicsof the structures; an output to a graphic system of nets, vector fields, isolines, field components, etc._ To verify the correctness of the operation of the program and its separate blocks, and also to estimate the accuracy, computation time, and the possiblity of calculating oscillations with multiple eigenvalues, problems which have an analytic solution wsre considered: oscillations in cylindrical and spherical resonators. Consider azimuthally homogeneous oscillations with the non-zero field components e,.e+,, h, in a spherical resonator. In a spherical set of coordinates the component h, has the form h,(r,O)=ll(kr)Pl(cosO), E---l, ii!,.... where it(z)are spherical Bessel functions, and P,(Z) are Legendre polynomials. The condition on the metal a(rb gives lower frequencies (in MHz), which are shown in Table 1, for a sphere of unit radius. The results of the numerical solution of this problem using the PIIUD-W system are shown in Table 2, from which we can obtain a representation of the convergence and accuracy of the method for a different number of nodes. In Table 2 the difference between the calculated and exact values of the frequency is denoted by S!, the linear dimension of the element is denoted by t&,andthe total time for calculating five frequencies on the computer with a speed of response of the order of 1Mflops is denoted by T (in seconds). The calculation time does not in practice depend on the initial approximation. In comparison withtheprogram from /l, 15/, where only the resonator problem is solved, the PRUD-W system provides now possibilities of obtaining the dispersion characteristics Fig.4 of periodic structures. This connected with the unconditional and rapid convergence of the space iteration method and with the possibility of obtaining immediately several (up to 101 eigenvectors. By specifying the phase shift and number of iterated vectors, we can obtain points immediately on several branches of the dispersion characteristics. The authors thank A.G. Abramov, S.YU. Ershova and T.D. Ryadov for their individual assistance and A.G. Sveshnikov for his interest. REFERENCES 1. WELLAND T., On the computation of resonances in cylindrically symmetric cavities - Nucl. Instrum. and Methods, 216, 329-346, 1983. 2. WEILAND T., On the numerical solution of Maxwell's equations and applications in the field of accelerator physics. Preprint DESY. Hamburg, 84-006, 1984. 3. KEIL E., Computer programs in accelerators physics. Preprint CERN. Geneva, 84-01, 1984. 4. MITTRA R. and LIE S., Analytic methods of wave quid theory, Moscow, Mir, 1974, 5. ANDREYEV V-G., Definition of the geometry of a structure with an alternating accelerating field in a x/Z--mode, Zh. tekhn. fiz., 41, 4, 788-796, 1971. 6. KEIL E., Diffraction radiation of charged rings movinginacorrugated cylindrical pipe Nucl. Instrum. and Methods, 100, 419-427, 1972. 7. ZOl'TERB. and BANE K., Transverse resonances of periodically widened cylindrical tubes with circular cross-section. - Preprint SLAC. Stanford, PEP-NOTE 308, 1979. 8. IL'INSKII A.S. and SLEPYAN G.YA., Oscillations and waves in electrodynamic systems with losses, Mir, Izd-vo MG, 1983. 9. NIKOL'SKII V.V., Variational methods for internal problems of electrodynamics,Moscow, Nauka, 1967. 10. RUTISHAUSER H., Computational aspects of F.L. Bauer's simultaneous iteration method. Numer. Math., 13, 4-13, 1969. 11. WILKINSON J.H., The algebraic eigenvalue problem, Moscow, Nauka, 1970. 12. WILKINSON J.H. and RAINSH S., Handbook of algorithms in Algal. Linear Algebra. tloscow, Mashinostroenie, 1976. 13. PARLET B., Symmetric problem of eigenvalues. Numerical methods. Moscow, Mir, 1983. 14. STRENG G. and FIX J., Theory of the Method of Finite Elements. Moscow, Mir, 1977. 15. HALBACH K. and HO&SINGER R.F., Superfish - a computer program for evaluation of RF cavities with cylindrical symmetry. - Particle Accelerator, 7, 213-222, 1976.

Translated

by H.Z.