Beam envelope calculations in general linear coupled lattices

The fundamental framework for the design and analysis of an uncoupled quadrupole lattice is the Courant-Snyder (CS) theory.¹ Many standard textbooks on accelerator and beam physics introduce the CS theory to begin discussions on linear uncoupled beam dynamics.^2–4 Interestingly, one of the recent trends in the beam physics community is to introduce coupling elements between two planes intentionally for sophisticated beam manipulations.^5–9 Several ideas have been proposed to generalize the original CS theory for the systematic beam dynamics description of such coupled systems with as much similarity as possible to the original CS theory.^6,10,11

In recent papers,^12,13 the dynamics of charged particles in general linear coupled lattices has been formulated in terms of a generalized CS theory, which extends the original CS theory to the fully coupled cases. The generalized CS theory can be applied to systems with quadrupoles, skew-quadrupoles, dipoles, and solenoidal components, as well as include the effects of torsion of the design orbit and variation of the beam energy. Furthermore, the generalized CS theory contains remarkably similar mathematical structure and physical meaning to their counterparts in the original CS theory. The 1D envelope equation is generalized to a matrix envelope equation in higher dimensions, and the phase advance is generalized to a symplectic rotation. Particularly in Ref. 13, the gauge group structure of the generalized CS theory was analyzed using the polar decomposition theorem. It has been shown that the gauge freedom in the matrix envelope equation for w can be bypassed by solving the symmetrized envelope equation for the generalized Twiss parameter β. This gauge adjustment makes the numerical algorithms of searching for matched solutions for β and $β′$ more efficient than for matched solutions for w and $w′$⁠. In the matrix envelope equation for w, in principle, there are an infinite number of matched solutions for w and $w′$⁠.

In this paper, we present a detailed analysis in Secs. II and III of the symmetrized envelope equation for β, and formulate the transfer matrix (M) and beam matrix (Σ) in terms of the generalized Twiss parameters α and β, which turn out to be remarkably consistent with the original CS theory. We also present beam envelope calculations in Sec. IV for several linear coupled systems, which indeed demonstrate the validity and usefulness of the generalized CS theory. Conclusions are summarized in Sec. V.

The general form of the time-dependent Hamiltonian in general linear focusing lattices is given by^12,13

where $z=(x,y,px,py)T$ are the transverse phase space coordinates. Here, $κ(s)$⁠, R(s), and $m−1(s)$ are time-dependent 2 × 2 matrices. Note that the matrices A_c, $κ(s)$⁠, and $m−1(s)$ are also symmetric. The superscript T denotes the transpose operation on a matrix, and $px(py)$ is the scaled canonical momentum variable conjugate to the transverse coordinate x(y) relative to the design orbit. The variable s plays the role of a time-like variable, and the prime denotes a derivative with respect to s.

If we apply the final results of the generalized CS theory described in Ref. 13 to the Hamiltonian in Eq. (1), we can construct the formal solution for the transverse beam dynamics in terms of a time-dependent linear map from the initial condition z₀ according to

where the transfer matrix M(s) is given by

Here,

The matrix u is the symmetrized version of the envelope matrix w used in our previous calculations,^14–17 and here it is equal to the square root of the generalized Twiss parameter β, i.e.,

Here, u(s) is symmetric and positive-definite according to the polar decomposition theorem. Because $β=wTw$ is also symmetric and positive-definite by definition, we note in Eq. (7) a well-known fact in fundamental linear algebra that a positive-definite matrix has a unique positive-definite square root.¹⁸ Therefore, u(s) can be uniquely determined by calculating the principal square root of $β(s)$⁠.

The symmetrized envelope equation for β is given by¹³ 

This is the main governing equation we need to solve numerically. Here, g and h are given by the lattice configuration as

The terms containing w and $w′$ can be expressed in terms of u, $u′$⁠, and D according to

Now the task is to express $u′$ and D in terms of β and $β′$⁠. In Ref. 13, two linear equations are derived as

We note that the above two linear equations are indeed special cases of a more general Lyapunov equation.¹⁹ While in Ref. 13 the inversion of a linear transformation L was formulated in terms of an eigenvalue problem, in the present work we use the general solution of the Lyapunov equation. The Lyapunov equation is a special class of linear matrix equations given by

where A, B, and C are assumed, for simplicity, to be square matrices of the same dimensions. Then the solution X can be shown to be

Here, I is the unit matrix, the symbol ⊗ represents the Kronecker product of two matrices, and the symbol “$vec$” represents an operator of a matrix which stacks the columns into a vector.¹⁹ Applying Eq. (16) in Eqs. (13) and (14), we can easily obtain expressions for $u′$ and D in terms of β and $β′$⁠. Now it is clear that Eq. (8) is a second-order nonlinear ordinary differential equation for the matrix elements of β. Making use of symbolic computational software such as MATHEMATICA,²⁰ we solve Eq. (8) numerically without having explicit expressions for the matrix elements.

To calculate the phase advance which is a 4D symplectic rotation, we need to solve the modified phase advance equation for P_u given by¹³ 

We choose the initial condition as $Pu0=I$ without loss of generality.

The transfer matrix in the original CS theory is decomposed into the elegant form² 

The third matrix on the right-hand side of Eq. (19) performs a canonical transformation of coordinates into normal forms, making the linear betatron motion a simple coordinate rotation given by the second matrix. The first matrix represents a back-transformation into the original coordinate system.

Since the governing equations in the generalized CS theory with gauge adjustment, i.e., Eqs. (8) and (17), seem rather complicated, at first glance it appears to be difficult to express Eqs. (4) and (6) in terms of Twiss parameters similar to Eq. (19). However, we note that if we define another Twiss parameter α as

we can express the transfer map into a surprisingly similar form to the original Courant-Snyder theory. We obtain

In deriving Eq. (21), we used the symplectic condition (⁠$SuJSuT=J$⁠, where J is the unit symplectic matrix) imposed on the canonical transformation matrix S_u, which can be expressed in terms of Twiss parameters as

This condition is trivial in the original CS theory. In the generalized CS theory, this symplectic condition is equivalent to calculating D properly. Indeed, Eq. (14) for D is derived from the symplectic condition. Therefore, once β and $β′$ are known from the envelope calculations, $α=−u(u′+Du−R)m$ can be calculated straightforwardly without violating the symplectic condition. This fact is also confirmed in our numerical calculations. It should be emphasized that the single-particle trajectory is determined in terms of the initial phase-space coordinates z₀ and lattice configurations, not from the form of the parametrization. Therefore, the transfer matrix itself is independent of the choice of β and $β′$⁠. Here, we have used the symbols “$β,α,M$” to represent physical quantities in both the original and generalized CS theories without causing any confusion. For example, β is a scalar quantity in the original CS theory, while it is interpreted as 2 × 2 matrix in the generalized CS theory. In this way, we can clearly display the correspondence between the two theories.

The beam matrix (second-order moments of the beam distribution) can be constructed by noting that the canonical transformation $z¯=PuSuz$ makes $z¯$ normalized coordinates in 4D phase space, and the normalized coordinates are invariant (⁠$z¯=const.$⁠) during the linear coupled motion.¹⁷ For the case of a Gaussian distribution, the beam matrix $Σ=〈zzT〉$ can be expressed as ¹⁷ 

where $⟨⋯⟩$ denotes the statistical average over the beam distribution. Here,

and ϵ₁ and ϵ₂ are two eigen-emittances. We note that $det(Σ)=det(ε)=ϵ1ϵ2≡ϵ⊥2$⁠. In terms of the Twiss parameters α and β, the beam matrix can also be expressed as

For the case of unequal emittances, the phase advance term does not cancel out. Hence the beam matrix depends not only on the Twiss parameters but also on the phase advance. When the two eigen-emittances are equal (⁠$ϵ1=ϵ2=ϵ⊥$ and $ϵ=ϵ⊥I$⁠), the beam matrix becomes

where γ is another generalized Twiss parameter defined by

This generalizes the familiar relation between the Twiss parameters α, β, and γ in matrix form. We note that the beam matrix is expressed in a remarkably similar form to the original CS theory.³ From Eqs. (3) and (23), we also confirm that the relation for transforming the initial beam matrix $Σ0$ is

In this section, we present numerical solutions of the beam envelopes for two lattice configurations using the generalized CS theory introduced Secs. II and III. The first case is a periodic lattice system with solenoidal focusing field,²¹ and the second case is the emittance transfer experiment (EMTEX) beam line at GSI, which is composed of solenoids, quadrupoles, and skew-quadrupoles.^7,22

We consider a Gaussian beam with negligible space-charge force propagating through a periodic focusing solenoidal field²¹ with axial periodicity length S = const. In this case, the symmetric matrices in the Hamiltonian (1) are $m−1=I$ and

The 2 × 2 matrix R is given by

where $Ω(s)=ebBz(s)/2γbmbβbc$ is the normalized Larmor frequency. Here, $Bz(s)=Bz(s+S)$ is the axial component of solenoidal field at $(x,y)=(0,0)$⁠, $βbc=const.$ is the average axial velocity, $γb=(1−βb2)−1/2$ is the relativistic mass factor, $eb(mb)$ is the charge (rest mass) of a beam particle, and c is the speed of light in vacuo.

When the applied focusing forces are periodic as in the periodic solenoidal system, it is desirable to determine the matched-beam envelopes. In general, the matched-beam envelopes provide minimum radial excursions in periodic lattice channels.²³ A general condition for the matched solution is²⁴ 

and it requires certain initial conditions. Here, for simplicity, we consider the case of equal eigen-emittances only. In this case, matched beam envelopes are found by imposing periodic boundary conditions when solving Eq. (8) as

In addition, we can compare our approach to the well-known method based on the Larmor-frame transformation.²¹ For negligible space-charge forces, the envelope equation for a circular cross-section beam with $σ=⟨x2⟩=⟨y2⟩$ is given by²¹ 

The matching conditions are

We choose a periodic step-function lattice with amplitude $Ω̂=const.$ and filling factor η illustrated in Fig. 1.

Plotted in Fig. 2 is the evolution of the matched beam envelope for one lattice period, $0≤s/S≤1$⁠. Here, the transverse dimension is normalized by $ϵ⊥S$⁠. It is clear in Fig. 2 that the numerical values calculated from the generalized CS theory are in excellent agreement with the numerical solution obtained from solving the envelope equation (33). This indicates that the formulation used in the generalized CS theory is valid, and also that the numerical routines for solving Eq. (8) are implemented correctly. The two envelope equations (8) and (33) are derived by two drastically different approaches. For example, in deriving Eq. (8) we do not involve the Larmor-frame transformation explicitly. While Eq. (33) cannot be easily extended to include more general coupled lattice elements such as skew-quadrupoles, Eq. (8) can be readily applied to systems with both skew-quadrupoles and solenoids.

For a single-pass system, the initial beam matrix $Σ0$ is given, and we only need to propagate it through the channel using Eq. (28). As mentioned earlier, the transfer matrix M itself is independent of the choices of β and $β′$⁠. Therefore, we can pick arbitrary initial values for β and $β′$ as long as they do not cause a singularity at s = 0.

We consider a beam line for the EMTEX facility at GSI.^7,22 The major components of the EMTEX section are a solenoid with stripper foil inside, a quadrupole triplet, and a skew quadrupole triplet. The parameters of the components included in this calculation are listed in Table I. The purpose of the EMTEX is to demonstrate transverse emittance partitioning. Under symplectic transformations, the eigen-emittances, determined from the beam matrix in canonical coordinates, are invariant. Therefore, to change the beam eigen-emittances, a non-symplectic transformation is required. In the EMTEX system, a charge stripper inside the solenoid performs the non-symplectic transformation. The charge stripping process changes the input beam rigidity $(Bρ)in$ to $(Bρ)out$ suddenly, and causes splitting of the initially equal eigen-emittances into two different values. The transverse inter-plane correlations are removed by a decoupling section following the solenoid. More detailed theoretical descriptions are given in Refs. 7 and 22.

We assume the initial beam matrix at the entrance of the solenoid is

The initial beam has equal horizontal (x-direction) and vertical (y-direction) rms emittances and no inter-plane correlation. Note that we adopt a canonical coordinate system $(x,y,px,py)$ instead of $(x,x′,y,y′)$ used in Refs. 7 and 22. In this numerical example, we consider a low-intensity beam of $D6+$ with energy 11.4 MeV/u. The beam is stripped to $3D2+$ in a 22 $μg/cm2$ carbon foil placed in the middle of the 0.2-m long solenoid. The initial beam parameters are $α0=0, β0=2.5$ mm/mrad, and $ϵ0=0.51$ mm mrad. During the stripping process, there is a finite amount of scattering $Δφ$⁠, which is estimated to be $Δφ=0.226$ mrad.²² The decoupling section described in Table I is set up for 1.0 T solenoidal field strength.

In Ref. 22, eigen-emittances and rms emittances along the EMTEX system are calculated by an analytical method based on the multiplication of a series of transfer matrices for each lattice element. Here, we apply the generalized CS theory, and calculate the eigen-emittances and rms emittances using a different approach. Since the envelope equation (8) (which is derived from the generalized CS theory) is valid for symplectic processes only, we need to divide the calculation domain into two regions: before stripping, and after stripping. First, we calculate the transfer matrix (⁠$Min$⁠) from the solenoid entrance to the stripping foil with initial beam rigidity $(Bρ)in$⁠. Then, we obtain the beam matrix at the stripping foil from $Σfoil=MinΣ0MinT$⁠. Finally, we apply another transfer matrix (⁠$Mout$⁠) from the foil with final beam rigidity $(Bρ)out$ to calculate the beam matrix along the decoupling section as $Σ(s)=Mout(s)Σ′foilMoutT(s)$⁠, where $Σ′foil$ includes the scattering effects on the angular spreads in $Σfoil$⁠. In this way, we can correctly calculate the propagation of the beam matrix across the non-symplectic element. It should be emphasized again that, different from periodic matching problems, there is no need to find initial values for β and $β′$ which yield the given initial beam matrix.

The eigen-emittances and rms emittances at the exit of the solenoid are plotted in Fig. 3 as functions of the solenoidal field strength. The eigen-emittances are calculated from

or equivalently,²² 

The results in Fig. 3 are obtained from the generalized CS theory, and they match almost perfectly with the results presented in Ref. 22. Partitioning of the eigen-emittances is clearly seen, and the initially uncoupled beam becomes strongly coupled. Once the beam passes through the EMTEX decoupling section, the inter-plane correlations are removed. Figure 4 shows the eigen-emittances and rms emittances at the exit of the skew quadrupole triplet. The decoupling processes are clearly shown, i.e., the eigen-emittances become equal to the rms emittances. The results in Fig. 4 reproduce the results presented in Ref. 22 almost exactly. It should be noted that the EMTEX decoupling section performs very well for a wide range of solenoidal field strengths.²⁵ While the solenoidal field strength is changed, the settings for the EMTEX beamline are kept constant at the values used to decouple the beam manipulated by a solenoidal field of 1.0 T. Figures 3 and 4 confirm the validity of the theoretical model and numerical procedure presented here, and also provide an independent cross-check for the emittance partitioning theory developed for the EMTEX system.

In this paper, we have formulated the transfer matrix (M) and beam matrix (Σ) in terms of the generalized Twiss parameters β and α. The recently derived symmetrized envelope equation for β, i.e., Eq. (8), has been solved numerically. The linear coupled beam dynamics is effectively described by the generalized Twiss parameters and a 4D symplectic rotation, which extends the elegant mathematical structure of the original CS theory. The validity of the formulation has been tested using two numerical examples. The numerical calculations based on the present formulation work very well for both periodic and single-pass systems. Furthermore, we note that the formulation can be applied to a non-symplectic lattice configuration by setting appropriate boundary conditions around the non-symplectic element. Since the formulation can include arbitrary linear coupled focusing elements as functions of s, it can be easily used to design and study a complex coupled lattice system for special beam manipulations, such as the EMTEX system.

This work was supported by the 2014 Research Fund (1.140075.01) of UNIST (Ulsan National Institute of Science and Technology). This work was also supported by the U.S. Department of Energy Grant No. DE-AC02-09CH11466.