Energy and RF Cavity Phase Symmetry Enforcement in Multi-turn ERL Models

In a multi-pass Energy Recovery Linac (ERL), each cavity must regain all energy expended from beam acceleration during beam deceleration, and the beam should achieve specific energy targets during each loop that returns it to the linac. For full energy recovery, and for every returning beam to meet loop energy requirements, we must specify and maintain the phase and voltage of cavity fields in addition to selecting adequate flight times. We find these parameters by a full scale numerical optimization program. If we impose symmetry in time and energy during acceleration and deceleration, fewer parameters are needed, simplifying the optimization. As an example, we present symmetric models of the Cornell BNL ERL Test Accelerator (CBETA) with solutions that satisfy the optimization targets of loop energy and zero cavity loading.


I. INTRODUCTION
An Energy Recovery Linac (ERL) reclaims energy from decelerating particle bunches to reduce the net power consumption of the accelerating radio frequency (RF) cavities. Return loops connect the exit end of the linac to the entrance, allowing the beam to recirculate through multiple accelerating passes of the linac. After a target highest energy is achieved, the beam decelerates through the same return loops. The energy of a decelerating beam is returned to the RF cavities, which reuse this energy to accelerate future particle beams [1]. The beam must satisfy two conditions: • Energy recovery. During deceleration, each cavity should regain the same amount of energy that it transferred to the beam during acceleration. If satisfied, the power load on each cavity from the beam will on average be zero.
• Reasonable energy targets. If the beam is intended for experimental applications, it must achieve the desired energy target. Additionally, the energy of the beam during each return loop must satisfy design criteria for the particular ERL construction.
Due to varying longitudinal velocity (v) of the beam throughout the ERL, it is challenging to find RF cavity phases that provide the desired acceleration during all linac passes. Adjustment of loop length alone may not guarantee the appropriate energy recovery in all cavities. If all cavity phases are adjusted to maximal energy gain for the highest energy particles, where v ≈ c, then low-energy beams will experience RF phases slipped away from maximal energy gain, and synchrotron radiation in the return loops can cause additional offsets in beam time of flights [1]. The individual phases and voltage settings of each RF cavity, as well as the time of flight through the return loops, are the parameters for optimizing cavity loads and the beam energies.
In the paper, we use a mixture of theory and numerical optimization algorithms and present a time-symmetric method of identifying loop lengths, RF voltage, and RF phase settings. This symmetry reduces the number of fit parameters required during optimization. We use a model of the Cornell BNL ERL Test Accelerator (CBETA), which has 4 physically distinct return loops and a linac that holds 6 evenly spaced accelerating cavities.
FIG. 1. Layout of the CBETA ERL [2]. A 6 MeV beam is generated in the injector, accelerates through a linac with six 7-cell cavities, circulates clockwise through the return loop, and recirculates to cover a total of 4 accelerating and 4 decelerating linac passes. The particle path concludes at the beam stop (top right). Splitter regions at either end of the linac independently control the flight times of beams with the target energies of the 4 loops.
CBETA brings a 6 MeV injected beam to 150 MeV over 4 accelerating passes and returns it to 6 MeV over 4 subsequent decelerating passes [2]. The 150 MeV beam, which travels around the 4 th ERL loop, can be used as a compact synchrotron radiation source or for internal target experiments in nuclear and elementary particle physics [3]. Afterward, the high-energy beam is recirculated through the full multi-turn path for deceleration.
As an example, if CBETA cavity phases are chosen purely for maximum acceleration in the first linac pass, the power load of a single transversely on-axis particle on each cavity ranges from 4 keV to 50 keV; particles with slight longitudinal or transverse offset from this ideal particle will incur even larger magnitudes of load per cavity. This load can be reduced by optimizing the phase set-tings. In this study, we develop a symmetric acceleration system with load and beam energy objectives. The system is implemented in CBETA simulations to illustrate the use of symmetry in enhancing the optimization process, and the resulting objectives are calculated in ERL models of increasing complexity.
We also generalize our approach to ERLs where the beam passes M times through a linac of N cavities. Each injected particle bunch results in M − 1 beams looping back to the linac at different points in time. In an ERL, the beam that has been accelerated by m ≤ M 2 linac passes has approximately the same energy as the beam that is yet to be decelerated m times; these two beams then traverse the same beam pipe, and one speaks of a M 2 -turn ERL.
If large energy aperture optics are used, all M − 1 returned beams can travel in one vacuum pipe. Such an optics configuration can be realized with a Fixed Field Alternating-gradient (FFA) design, as in CBETA [2]. Assuming the beam energy stays constant during the return loop, then the energy changes occur a total of M ·N times over all cavity encounters.
If the profile of energy gain and loss is symmetric through the entire ERL, then we call this a symmetric ERL. With this symmetry, the energy that the beam gains when it traverses a cavity for the i th time is the same as the energy that the beam loses during the (M N − i + 1) th time.
Using these notations, CBETA is a 4-turn ERL with N = 6 and M = 8. Because of its FFA optics, CBETA only has one vacuum pipe for all 7 returned beams. Additionally, by adjusting parameters appropriately, it can be operated as a symmetric ERL.

II. OPTIMIZATION SYSTEM PARAMETERS
Suppose an M 2 -turn ERL has N cavities and M linac passes. To indicate individual loops and cavities, we use indices m and n such that 1 ≤ m ≤ M and 1 ≤ n ≤ N .

A. Objectives
During ERL operation, beam timing must ensure that an appropriate amount of energy is transferred during each and every cavity passage. Optimal energy recovery requires there to be a minimal load on all RF cavities: where ∆E mn is the energy gain in cavity n during pass m. When optimizing an ERL, we use E load,n as N objective functions that are to be minimized. For the beam to follow the design orbit, the beam energy E loop,m at the end of the m th linac pass must be close to the design energy, E des,m , to be compatible with the magnet settings, where E loop,0 is the beam energy at initial injection. An optimization then minimizes the objective functions |E loop,m − E des,m | for 1 ≤ m ≤ (M − 1). In a perfect ERL where E load,n = 0 for all n, Eq. (1) leads to E loop,M = E loop,0 . Therefore, the objective function for m = M is not used, and we have N + M − 1 objectives.
A symmetric ERL has E des,m = E des,M −m . If the operation of the cavities guarantees a symmetric energy profile, then also E loop,m = E loop,M −m during each optimization step, and the beam loading is symmetric as well: E load,n = E load,N −n+1 . If N is odd, then the central cavity in the linac already has zero load when ERL symmetry exists, and the central load ceases to be a useful objective. Using the Gauss bracket to designate the floor of a real number, we then obtain only [ N 2 ]+ M 2 objectives. One can further reduce the number of objectives by not requiring specific target energies in the intermediate loops, and only specifying the highest energy of loop M 2 . The magnet settings would then be adjusted accordingly to the resulting energies of the intermediate loops. This reduces to [ N 2 ]+1 objective functions. For CBETA, we therefore generally have 13 objectives. As a symmetric ERL, we have 7 objectives. When intermediate energies are not fixed, we only have to meet 4 objectives.

B. Degrees of Freedom
To identify degrees of freedom that will optimize load and loop energy objectives, we need independent parameters that affect ∆E mn . In general, we can control the following parameters that affect cavity energy gain: • Initial RF phase of the n th cavity, φ 0,n .
• Voltage of the n th cavity, V n .
• Travel time through the m th loop, t loop,m .
These generally provide 2N + M − 1 degrees of freedom.
In a symmetric ERL, the last [ N 2 ] cavities are operated symmetrically to the first [ N 2 ], and therefore do not provide phase or voltage degrees of freedom. If N is odd, then the phase of the central cavity must also be chosen such that it operates symmetrically to itself; however, the voltage of this cavity is still free. Together, the [ N +1 2 ] voltages and [ N 2 ] phases yield N degrees of freedom from the cavities. The energy symmetry also makes the (M − m) th time of flight equal to the m th one. There are then N + M 2 degrees of freedom. In the case of CBETA, we generally find 19 degrees of freedom. With fully symmetric operation, we have 10 degrees of freedom.

C. Efficient Optimization
In general, we can make the optimization more efficient by having a small number of objectives and equally few degrees of freedom. In each ERL system, there are more degrees of freedom than objectives. To speed up the optimization process, one can decrease the number of degrees of freedom. For example, one can set all cavity voltages to the same value. In a general ERL, there are N more degrees of freedom than constraints; we can therefore make the N voltages constant to achieve an optimization system with an equal number of objectives and degrees of freedom.
A During optimization, a solution generally cannot be found if there are fewer degrees of freedom than objective functions. For example, the degrees of freedom for the symmetric CBETA system should not be minimized with both constant V n and constant t loop,m = M 2 , because this will result in using only 4 degrees of freedom to fulfill 7 objectives (Table I, bottom row).

III. ERL SYMMETRY CONDITIONS
A symmetric ERL has a significantly reduced number of objective functions and degrees of freedom compared to a non-symmetric ERL, and it can therefore be optimized much more easily. This is true for accelerator simulations as well as for experimentally finding the desired accelerator settings. In this section, we therefore discuss operation conditions that create a symmetric ERL. This is achieved most easily if the cavities are arranged with mirror symmetry about the center of the linac cryomodule. This means that the distance from the n th cavity to the next is the same as the distance between the N − n th cavity and its next neighbor. There will also be cavity phases where the standing waves in the n th cavity have a field that is the mirror image of the (N − n + 1) th cavity field. Mirror symmetry in the fields always exists if the geometry of the n th cavity is the mirror image of that of the (N − n + 1) th one.
In CBETA, the cavities are designed identically, but they are installed such that all are facing the same direction: an input coupler is located on the downstream end of each cavity. However, because the cavity cells are constructed in a symmetric way, and the field from the input coupler is small compared to the main cavity field, these cavities are symmetric to a sufficient degree.
Mirror symmetry in the cavities, as well as over the distance between cavities, can result in a symmetric ERL if we choose the particle timing and cavity phases such that following a particle forward from injection, or tracking it backward in time from the beam stop, both result in exactly the same electric fields seen by the particle. If these conditions are met, the acceleration and deceleration profiles of the forward and backward-traveling particle are then identical, and we have a symmetric ERL, i.e. the beam decelerates with the same energy steps as during acceleration.
To produce a symmetric ERL, we need a way to set cavity phases such that accelerating and decelerating cavities can be symmetric. In the following discussion, we first show that symmetry can be established in a linear sequence of two cavities. We then transfer this symmetry to the closest equivalent ERL: a 2-pass, 1-cavity system, where the same cavity is traversed twice. The ERL is then expanded to an arbitrary even number of M passes and 1 cavity. Finally, we combine the findings of the linear and 1-cavity systems to arrive at the symmetry requirements for a general M -pass, N -cavity ERL.

A. Linear Sequence: 2 Cavities
Imagine a straight linac with two independent cavities (A and B) of frequency ω, arranged in a mirror symmetric way about the center between the cavities. We can represent the two cavities as mapping functions C {E, φ}, The goal is to find an input phase when the particle enters the second cavity, φ in,B , such that the first cavity adds as much energy as the second cavity subtracts. We define the cavity output phase, φ out,A , as the input phase plus the particle's travel time through the entire cavity, The input phase of the second cavity, φ in,B , must be chosen such that the output phase, φ out,B , produces a field that is equal and opposite to the field when the particle enters cavity A, while the time derivatives are the same. The field of a standing wave in cavity A depends on longitudinal position s, time t, and on the input phase. Between start position s = 0 and end s = L of the first cavity, the field observed by the particle is, where the particle enters the cavity field at the input time t = t in,A . The cavity input phase φ in,A is independent of particle entrance energy, but it can be expressed in an energy-dependent form for particle speed v ≤ c, where the relative input phase, φ v in,A = 0, is the parameter commonly controlled during typical ERL operation. The phase offsetφ v A enables the particle with speed v to experience maximum acceleration at an on-crest phase of φ v in,A = 0. In general,φ v A =φ v B because these offsets depend on the energy of the particle upon cavity entrance; however, an ultra-relativistic particle will result in identicalφ c for symmetric cavities A and B.
The spatial RF field dependence is given by E A0 (s), and by convention it is chosen to start with a positive value in the first cell. For symmetric cavities with an odd number of cells, the spatial dependence is a symmetric function about the center; for an even number of cells, it is an anti-symmetric function. This means that E B0 (L − s) = ±E A0 (s). The sign (+) is for odd and (−) is for even numbers of cells per cavity.
The exit of cavity B is placed a distance L s after the entrance of cavity A, and E describes the net field that the design particle experiences from the two-cavity sequence, including A and B contributions. The design particle travels through these cavities so that it reaches location s at the time t(s). In cavity A, therefore E (s) = E A (s, t(s)). At position (L s − s) within cavity B, the particle should experience an equal but opposite field as during acceleration in cavity A at location s, where Eq. (7b) follows from Eq. (7a) if the specific cavities are considered. The time spent from position s to the end of cavity A must be the same as the time spent from the start of cavity B to a distance s before its end: where T A is the full duration that the particle spends in cavity A, as previously defined. We then equate the mirror symmetric field, E B (L − s, t(L s − s)), with the cavity A field in Eq. (5) using the field symmetry relation from Eq. (7b) and time from Eq. (8), Eq. (9) leads to the cavity B phase input, φ in,B , for either odd or even-cell cavities, using φ out,A = φ in,A + ωT A , If the input phase of B is chosen according to this relation, then traveling backward in cavity B will automatically give the same time-dependent fields as traveling forward in cavity A.
To verify in a simple example that acceleration is indeed reversed for an ultra-relativistic particle in symmetrically phased cavities, consider cavities A and B as pillboxes operating in the fundamental mode with β = 1, where one pillbox is traversed in half an oscillation by a particle with v = c. The energy gain of this ultrarelativistic particle is given by, where q is the charge and V A = V B is cavity voltage; since the particle is ultra-relativistic, the qV A pre-factor can be taken as a known, identical constant in both cavities [4]. If A and B are each sets of neighboring pillboxes in a multi-cell π-mode cavity, then the ultra-relativistic particle crosses cavity A with a time of flight ωT A = n cell π, where n cell is the odd or even number of cells per cavity. To find the energy change from B, apply the Eq. (10) symmetry conditions to the input phase of B, This verifies the Eq. (10) symmetry conditions for odd and for even cell numbers: cavity B subtracts as much energy as cavity A adds.
We now search for a way to establish similar symmetry in an ERL. Consider an ERL with the minimum number of cavities and linac passes, N = 1 and M = 2. A particle in this small ERL will encounter this cavity twice: first when accelerating (encounter 1), and secondly when decelerating (encounter 2). Unlike in the case of two successive cavities, the input phase for the second encounter, φ in,2 , depends on the original cavity phase choice, where t pair is the amount of time spent between leaving the first encounter and entering the second. If the cavity is symmetrically built and the proper t pair is found to satisfy the linear sequence symmetry condition from Eq. (10), the 1-cavity ERL will become symmetric. Substitute this condition into Eq. (12), setting φ in,B = φ in,2 , to determine the necessary t pair for symmetry, (13)

C. ERL: 1 Cavity, M Passes
We now extend the symmetry to an ERL with one cavity, but M ≥ 2 passes. Let m represent the pass index, where 1 ≤ m ≤ M . The m th encounter of the cavity must hold a symmetric phase relation with the (M − m + 1) th encounter, in the manner of Eq. (10).
In the form of Eq. (13), t pair,m designates the time spent from the exit of the m th cavity encounter to the entrance of the (M + 1 − m) th encounter. For each m, this time between pairs must be correctly related to the phase of the m th encounter by Eq. (13). This is done by adjusting the time spent in the m th return loop, t loop,m , between the exit of the m th encounter and the entrance of the (m + 1) th encounter.
The time between cavity pairs, t pair,m , includes the times T j spent within the cavity of encounter j, for j = (m+1) through (M −m), as well as the in-between t loop,k times of flight, where k = m through (M −m). Therefore, where the formulation takes into account that T k = T (M −k+1) and t loop,k = t loop,(M −k) due to symmetry. Conversely, the intermediate loop times become, for m < M 2 , and the highest energy loop has t loop, M With the full t pair,m expansion from Eq. (13), this leads to the loop times, This can also be expressed in relative phases, φ v in,j = φ in,j −φ v j , instead of absolute the φ in,j notation, for any encounter j; then, in the form of Eq. (6),φ v j relates the input phase to the RF field profile required for maximum acceleration at φ in,j = 0. If the particle is ultrarelativistic, then all phase offsets are identicallyφ v j =φ c , and t loop,m from Eq. (16) can be written explicitly in relative phases, This example shows that the phase difference, as one would expect, simply accumulates during the travel time if the reference phasesφ v j do not change between cavity encounters.

D. ERL: N Cavities, M Passes
Consider a large ERL with N cavities and M passes. As in the 1-cavity ERL, the input phases φ in,mn of all cavity encounters after the first pass are dependent on the initial encounter phase, φ in,1n , where again 1 ≤ m ≤ M and 1 ≤ n ≤ N . Drift pipes between consecutive cavities increase particle times by t drift,mn between the n th cavity and its subsequent neighbor, while the return loops have time of flight t loop,m . The particle is injected into the ERL at a time t = 0, and it enters cavity n on pass m at the input time t in,mn .
Let t total represent the full time that a beam would take to travel through the symmetric ERL. For brevity, n = N − n + 1 and m = M − m will designate the cavity/pass indices for the cavity that is symmetric to the one with n and m, where we let n , m refer to the earlier encounter and n, m to the later encounter, i.e. m > m . As in the example of two successive cavities, the voltage of n and n must satisfy V n = V n to permit mirrored cavity interactions during symmetric acceleration and deceleration. If N is odd, then n = n for the central cavity, but this case can be treated identically to all other paired cavities.

Pass 2 Deceleration
Symmetry Axis ϕ in,12 t=t in,12 ϕ out,12 t=t out,12 FIG. 2. Cavity encounter sequence for a 2-pass ERL with N = 2 physical cavities, denoted C1 and C2, and cavity phase inputs φin,mn and outputs φout,mn associated with corresponding input/output times t. The total time spent in the ERL, t total , is represented by the length of the horizontal line. This schematic can be extended to general multi-pass ERLs by inserting additional cavities or return loops on both sides of the dotted symmetry axis.
If an input phase φ in,mn is known, we can use it to define the effective initial cavity phase. Let φ 0,n be the total input phase of a particle that enters cavity n at time t = 0, where the beam exits the cavity at time t out,mn = t in,mn + T mn . The phase of the n th cavity must satisfy Eq. (10) to reverse the primed cavity acceleration. We then substitute φ in,mn and φ out,m n from Eq. (18) into the symmetry conditions in Eq. (10) to solve for the unknown φ 0,n in terms of time and known, primed quantities, When the ERL is symmetric, as we seek to achieve, any point in time from the Fig. 2 sequence can be described by either counting in equal steps forward from the beginning, at t = 0, or backwards from the end, when t = t total . Cavity entrance times follow this principle, If t total is taken as a known parameter, we can relate the primed and unprimed initial phases. Substituting Eq. (20) into Eq. (19) for even or odd cells, the φ 0,n phase constraints for ERL symmetry become, To ensure that dependent phases take the proper values in passes m > 1, the parameter t total must accurately describe the total time that a beam will spend in the ERL. Since the highest energy loop crosses the symmetry axis of the unwrapped ERL, it is easiest to get a correct t total by setting the length the M 2 th loop to achieve an appropriate time of flight, where t out, M 2 N is the full duration that the beam spends accelerating, from particle injection to exiting the last cavity of the highest-energy pass.
When adjusting the cavities of a physical ERL, it is typical to only have control over the velocity-dependent φ v in,1n relative input phases during the first pass. Symmetry conditions formatted for direct ERL control would then require rearranging Eq. (21) to solve for the firstpass φ v in,1n relative phases, [even] (23) If the particle is ultra-relativistic, then paired cavities n and n have the same phase offsetφ c n in all passes, and Eq. (23) simplifies, For complete ERL symmetry, all dependent φ 0,n initial phases with n ≥ n must be constrained according to Eq. (21), which may be implemented as either Eq. (23) or Eq. (24) depending on particle speed. The M 2 th flight time must also follow Eq. (22). Intermediate return loops may then be set at any reasonable time of flight without disrupting this symmetry.

IV. CAVITY MODELS FOR ERL USE
The symmetry conditions and objectives are tested in ERL simulations that use one of the following cavity representations to calculate particle energy and time of flight interactions. Each model is an approximation the real CBETA cavities, which have 7 cells with an elliptical geometry [2].
• Ultra-relativistic Cavities (UR). Particles are treated with ultra-relativistic flight time and energy change formulas inside multi-cell pillbox cavities.
• Finite Time-tracked Cavities (FT). Cavities are single-cell pillboxes where non-ultrarelativistic particle behavior is considered.
• Runge Kutta Cavities (RK). Non-ultrarelativistic particle time and energy effects are calculated via integration through multi-cell pillbox cavity fields.
All models track particle time and energy within cavities with ideal pillbox cell shapes. The first three variants are tested in custom Mathematica programs. The latter, which is used as the accuracy reference for the other three, is modeled in the Bmad accelerator toolkit [5]. Each model is compatible with the symmetric and reduced symmetric ERL objective systems from Table I.

A. Thin Lens (TL) Cavity Model
As a first-order approximation, RF cavities can be modeled as a delta-function acceleration over infinitesimal time and distance. Each cavity imparts an instantaneous energy kick consistent with an on-axis beam in an infinitely thin, pillbox-shaped cavity, where the relative input phase, φ c in , is used in the same form as Eq. (11) due to the negligible cavity transit time. The voltage V describes the maximum possible energy gain of a q-charged particle, at the relative input phase φ c in = 0. Thin lenses have the position and entrance time of the midpoint of a typical finite-length cavity. Since these cavities have zero length, and effectively zero cells, symmetry in time must be determined using the phase constraint from Eq. (10) for even-cell cavities. Drift pipes affect the particle time according to particle velocity (t drift = L drift v(E) ). The pipes on either side of a thin lens are extended by a distance that compensates for the missing length of the thin cavity, thereby keeping the same overall linac dimensions as a more realistic system.

B. Ultra-relativistic (UR) Finite Cavity Model
The thin lens model neglects the length of a physical cavity. For a better approximation of ERL interaction, we modify the TL equation to model beam energy and time behavior as ultra-relativistic within a cavity of nonzero length. Outside of a cavity, the beam still has typical energy-dependent velocity. For example, if a particle of charge q travels through a 7-cell π-mode cavity (n cell = 7) of length L = (n cell πc)/ω and voltage V , then time and energy are modeled as, This energy change and time of flight across the cavity are only accurate for a particle with speed c. The general energy gain of any traveling particle with original energy E in would follow a more complex relation, and this relationship is further explored the FT model. The UR model assumes that the electron beam travels at the speed of light, which is not true for a typical MeVorder ERL injection energy. We take this into account by approximating that the particle spends the first half of the cavity length with the initial velocity and half the distance with its final velocity: with velocity calculated as, we then find the momentum, Inserting T FT into Eq. (29) leads to a transcendental equation for p out that can be solved numerically, from which the final energy is obtained by E out = c (mc) 2 + p 2 out . The energy difference is then, where E in is the original particle energy when it enters the cavity.

D. Runge Kutta (RK) Cavities
In the Bmad accelerator toolkit [5], cavities can be tracked with the Runge Kutta 4 algorithm. Particle energy and position are found via integration through the electric and magnetic fields of the cavity. If the particle travels through the exact center of the cavity (transversely on-axis), then it will only experience a longitudinal electric field, E s . If the cavity is an ideal pillbox, and the Bmad auto-phase calculation is deactivated (direct φ in control), then Bmad models the longitudinal s component of standing wave fields as, where the particle encounters voltage V over the pillbox length L = πc ω for RF frequency ω, the speed of light is c, and the wave vector is k = c ω [5]. Note that for a pillbox run in the fundamental harmonic, k = 0. The field in Eq. (31) is effectively the first harmonic of a pillbox cavity. Since first harmonic fields are anti-symmetric about the center of a cell, but CBETA cavities operated in the fundamental harmonic should be symmetric about the cell center, the typical first harmonic spatial and time dependencies are shifted by π 2 to create a center-symmetric field pattern. This modified first harmonic field models the physical field in elliptical cavities slightly better than the fundamental pillbox mode would.
Cavities with multiple cells are represented as series of consecutive pillboxes stacked end-to-end. The particle enters the first cell with input phase φ in ; the output phase of the first cell becomes the input phase for the second cell, the second output feeds into the third input, and so forth. This setup approximates the effect of a π-mode cavity, although it does not perfectly account for nonultrarelativistic time dependence between different cells.
When modeling symmetry in CBETA with RK cavities, an on-axis particle is tracked through a 7-cell RK pillbox using a 1D Runge Kutta algorithm.

E. Cavity Model Accuracy
Physical CBETA cavities have 7 cells with an elliptical geometry. For the purposes of this study, we will approximate these as 7-cell pillbox stacks. To determine how well the mathematical TL, UR, and FT approximations model a real 7-cell pillbox system, we compare them to the most rigorous model: a RK cavity. The RK model is the most realistic in this regard because it directly integrates the particle through the electric fields, while the other models use approximations to determine the travel time and particle energy profiles.
For this comparison, we introduce the FT7 model, which represents a particle moving through 7 identical FT pillbox cavity cells. This is an extension of the FT model, which in its original form describes a single pillbox. The FT7 behavior is included as a demonstration of the difference in accuracy between 1-cell or 7-cell time and energy effects, but FT7 is not used in any subsequent full ERL models.
In Fig. 3, an on-axis particle with an energy of 6 MeV at the cavity entrance is tracked through TL, UR, FT, FT7, and RK models. Cavities are designed with the CBETA frequency ω = (2π · 1.3 GHz), and one cell has length πc ω . Cavities with an active region shorter than 7 cells, such as the zero-length TL model or the 1-cell FT model, are extended to a comparable 7-cell length by adding symmetric drift pipe extensions to either side of the active cavity. Time of flight is then measured from s = 0 to s = L, where L = 7πc ω is the length of a 7cell cavity. In all models, phase is defined as the φ in at the beginning of the cavity itself, not the drift pipe extensions.
The TL and UR models have an identical, substantial departure from RK energy gain, with errors on the order of 100 keV. The UR time of flight is also the least accurate, due to its assumption of particle v = c within the entire cavity; the TL and FT have near-identical time of flight behavior because the extended drift pipes on either side of the TL cavity create a similar time-average effect as the FT time-average calculation in Eq. (27).
Of the models considered, FT7 is most accurate to the RK results, and therefore the closest to expected physical cavity performance. If accuracy to physical systems is the primary concern, then the FT7 should be used; however, if simulation speed is also important, we should consider a trade-off between the desired criteria. The FT7 model requires the simulator to numerically solve 7 sets of Eq. (27) and Eq. (30), which results in a more time-consumptive calculation than the original FT model. This is not ideal, particularly when a numerical optimization may require thousands of iterations to converge on a solution with the desired objective precision.
For our purposes, the differences in time or energy accuracy between the FT and FT7 models are negligible (Fig. 3) compared to the increase in simulation speed. In the interests of running a large number of simulations within a reasonable program run-time, the full ERL models in this study consider only the TL, UR, FT, and RK cavities.   6.0000 6.0000 6.0000 6.0000 Optimized values of objective functions, defined in Eq. (1) and Eq. (2) as the differences between modeled load or energy and respective targets, are provided for each model in Table II. The degrees of freedom and fixed input values are also present.
The optimal φ 0 and t total settings for UR, FT, and RK models correspond well, but the optimized TL phase solutions deviate by around 10-15 • by model. This may be because φ in for TL cavities is measured at the start of the thin lens itself, which is located in the center of a nonzero-length UR, FT, or RK cavity. The slight difference in intermediate t loop return time of flights may also be related to this phenomenon, as the values were chosen for each model individually to achieve pre-optimized beam energies close to the design energies.
The three Mathematica models arrive at approximately machine precision on cavity load and central loop energy objectives, as defined in Eq. (1) and Eq. (2). For comparison, the typical beam energy falls on the order of 10-100 MeV, and the objective offsets are on the order of 10-100 µeV.
The 2-3 order of magnitude lower precision on the RK solution is a consequence of the Bmad numerical optimization algorithm used; the convergence of the lmdif Bmad optimizer on a solution depend highly on the initial parameter values, objective weights, and step size chosen for optimization. The larger values of the load 4-6 objectives in the RK model, which were not explicitly used in optimization, appear because the RK model setup only achieves about 0.1 eV precision between symmetric accelerating and decelerating energies, as opposed to the µeV-scale energy symmetry established in the three Mathematica models. Nevertheless, in this solution, the RK objective function most poorly satisfied (peak energy, ∆E loop,4 ) is less than 1 eV away from the 150 MeV target.
With the introduction of ERL symmetry in the four models of increasingly complex cavities, the necessary optimization system is reduced to a 4-by-4 set of objective functions and degrees of freedom. Although only half of the cavity loads were used in optimization, the solutions from Table II indicate that symmetry does result in correspondence between the optimized loads 1-3 and the non-optimized values of cavities 4-6. These solutions therefore show the practical use of ERL symmetry to satisfy objective functions.

A. Solution Sensitivity
Once optimal solutions are identified, we introduce errors to the input parameters to determine the accuracy required for successful ERL operation. The inputs to be varied include all N phases, N voltages, and M 2 loop lengths in the system; output functions are the objectives of maximum beam energy and all N cavity loads. For small input permutations about the Table II solutions, the objective function response is approximately linear; we can therefore speak of an approximately constant slope (sensitivity) of each objective to a single degree of freedom.
Objectives must be kept within a certain tolerable range around the ideal values (zero load, 150 MeV peak energy) for the ERL to operate in an acceptable manner. Assume that all inputs except one, denoted j, are set to the optimized values from Table II. We denote each objective function as f (j) for a particular objective f , such that f (j) = E load,n (j) or f (j) = ∆E loop, M 2 (j). The input j is therefore one of 2N + M 2 total input parameters, and the function f (j) is one of the N + 1 objectives in consideration.
The objective tolerance, f 0 , can be divided by that same objective's sensitivity to error in j, d dj f (j), to find the maximum error ∆j f that the input j can have before f (j) exceeds tolerance. The smallest value in the set of {∆j f }, for all functions f (j), therefore satisfies all f 0 tolerances.
If all inputs are at optimized settings except one, designated j, then we can vary input j from its optimized solution until the magnitude of objective function f (j) reaches its tolerance limit. The CBETA f 0 tolerances are chosen as a maximum load of 50 keV per cavity (2 kW for a 40 mA injection current), and 150 keV offset from peak beam energy. We then use Eq. (32) to determine the ∆j of each input when only that input is varied from the ideal solution. In Table III, the TL phase and loop length inputs show the lowest sensitivity, while the RK model phase and length inputs are the most sensitive to individual error sources. All models have approximately the same voltage sensitivity. The sensitivity results are useful to establish upper bounds for individual input error. However, in reality, the combined effect from multiple error sources must be considered. Further work is required to study the combined sensitivity. Suppose that we only have control over the input designated j. Furthermore, suppose the expected error on every other input i, for i = j, has a known standard deviation, σ i0 . Using statistical error propagation formulas, we can calculate an allowed error range σ j f for the special input j, such that the objective f (|j| ≤ σ j f ) is still within tolerance,

(33)
Error ranges for all objectives f (j) can be combined into the set σ j f , which has as many elements as the number of system objectives. The smallest σ j f is the range of j required to satisfy the most restrictive objective; therefore the range required to satisfy all objectives is σ j = min({σ j f }).
Using the same tolerance bounds as for the singlevariable errors, we set σ f 0 = 50 keV for load, and σ f 0 = 150 keV for the maximum beam energy. Expected σ j0 fluctuations for phase, voltage, and return loop length inputs are 0.1 • , 600 eV, and 1 3 mm, respectively. A calculation of combined tolerances using Eq. (33) gives the Table IV error TABLE IV. Error tolerance ranges for cavity phase, voltage, and return loop length, as calculated using Eq. (33) with expected σj0 fluctuations on all input parameters. For a system that can feasibly achieve and maintain operation within predefined objective bounds, the input tolerance should be larger than expected σj0. Table IV represent the stability required from the cavity field regulation system and the resolution required from the path length adjustment system. If errors on the order of σ i0 exist on all inputs but the one denoted j, Table IV indicates that input j must be regulated about twice as strictly as the others to prevent objective values from exceeding the σ f 0 tolerance bounds.

VI. LONGITUDINAL BEAM DYNAMICS
The modeled solutions and sensitivities thus far have only considered a single-particle, on-axis beam with zero time or energy spread. A more realistic particle distribution is expected to have non-zero dimensions in phase space (x, y, z, p x , p y , p z ), where all coordinates are relative to the position of the ideal single-particle scenario. Since the models only track the particle in the longitudinal phase space, consider a transversely on-axis beam (< x, y, p x , p y >= 0) with finite spread in injection time and energy (momentum). All longitudinal phase space coordinates are defined with respect to an ideal particle with coordinates z = 0, δ = pz−p0 p0 = 0. The injected beam follows a Gaussian distribution in time and energy, with σ time = 4 ps (σ z ≈ 1.1 mm) and σ δ = 5 · 10 −4 . In the desired scenario, an injected beam should achieve minimum δ spread during its highestenergy pass. This will allow better control over the individual particle energies after acceleration if a fraction of the beam is siphoned off for experimental purposes. Furthermore, to preserve the acceleration-deceleration symmetry found in the optimized ideal particle solutions, it is preferred if the beam has equivalent energy spread in the longitudinal phase space profiles at the injector and stop.
When a Gaussian beam is run through an ERL with settings optimized for the ideal particle, the longitudinal phase space profile becomes distinct from the original distribution after a complete, 8-pass ERL circulation.   4. Time evolution of 1000-particle Gaussian beam in longitudinal phase space: an initial Gaussian (both plots, overlaid), the beam after complete acceleration (left), and at beam stop (right). The beam is in an FT model ERL with M = 8, N = 6, and the ERL is set with ideal particle solutions (Table II). The ideal particle is defined at coordinates (z, δ) = (0, 0). Vertical and horizontal lines represent one standard deviation of the stated output from the mean z or δ coordinate.
Note that even in a symmetric ERL, the phase space distributions of a non-ideal beam before injection and after beam stop are not necessarily symmetric. The high energy spread at beam stop may be reduced if more objectives and degrees of freedom are added to the optimization system. A new objective could limit the energy spread of the final beam to a certain tolerable range, but this would also require new degrees of freedom to ensure a useful optimization system. For example, the voltages of paired cavities could be independently varied instead of set to a common value.
An easier solution to the phase space asymmetry involves adjustment of the injected particle distribution. Since the ERL converts a flat Gaussian z, δ distribution into a diagonal one, it is reasonable to suspect that, if we inject a beam with initial tilt in the opposite direction, a symmetric mirror image of the phase space distribution may emerge after the ERL (e.g. Fig. 5, right). This beam would then have equal energy spreads in injection and stop. Such an injection pattern may also decrease the energy spread during the M = 4 return loop, which will allow better control of the beam for experimental purposes.
The amount of tilt can be adjusted by varying the phase of the final cavity in the injector section. A proper choice of tilted input beam, such as the example in Fig.  6, allows equal initial and final σ δ .
To identify the correct input distribution for a symmetric output energy spread, consider the time and energy effects across the complete ERL, from injection to beam stop, as a longitudinal transfer matrix, Elements can be evaluated numerically by injecting test particles at known small (z inj , δ inj ) offsets from the ideal particle, which is defined at coordinates (0, 0), and measuring the final (z stop , δ stop ) values at the end of the ERL. We additionally construct an equivalent [R ij ] matrix, which describes the longitudinal phase space mapping of a particle that travels backward through the ERL, from injector to beam stop, Due to the symmetry between ERL acceleration and deceleration, the mapping of (z inj , δ inj) to (z stop , δ stop) must be identical in both Eq. (34) and Eq. (35) in order for particles with offsets (z inj , δ inj) = (0, 0) at injection to regain offsets of the same magnitude at the beam stop. Therefore, R 11 = M 22 , R 12 = M 12 , R 21 = M 21 , and R 22 = M 11 .
If the initial particle has only a small offset from the ideal case, then the ERL mapping can be considered linear, and the two coordinates are linearly related (δ = Az). For the initial and final beam distributions to have equal and opposite tilt, it is required that, An injected beam with small σ z offset must then have the net linear slope, With this slope A, the beam will exit the ERL with the same σ δ energy spread and an opposite tilt angle as it had at injection. In Table V, a particle with z = 10 −7 m or δ = 10 −7 offset is used to calculate slopes with Eq. (37). These slopes are used to add a vertical offset to all particles in a Gaussian distribution of (z, δ) particle coordinates (Fig. 5); we then report the tilt amount as the angle of a linear fit of the beam from the positive z-axis.  For beams with large σ z , the linear matrix transform no longer accurately describes the behavior of particles at the edges of the beam. With this type of beam, the proper tilt pattern is determined by injecting multiple test beams with different injector cavity phases, φ in .
A scan for the proper initial beam is conducted by injecting multiple test beams of different longitudinal profiles through the ERL. First, we simulate the pre-ERL injector cryomodule as a single accelerating UR-model cavity with a voltage of 1.5 MV. A beam with Gaussian energy and time distributions passes through this injector module. This beam is distributed about a pre-injection particle energy and time, E 0 and t 0 , such that after tilting, the central particle will have ideal particle injection properties, E inj = 6 MeV and t inj = 0 s. Pre-injection parameters E 0 and t 0 are determined by the reverse process of calculating which input energy and time will result in ideal injection parameters after accelerating through the injector module. The tilted beam is sent into the ERL, and its output σ δ is measured.
We then repeat the process with a new test beam by varying the injector phase in steps of, for instance, π 100 , while using the same pre-injector Gaussian distribution. Once the full range of phases is covered, from 0 to 2π, we identify the injected beam that yields the most similar σ δ energy spread in the injected and output distributions, where σ δ,stop σ δ,inj → 1. The ideal tilt would yield equal energy spread at injection and beam stop, but practical beam results are limited by the resolution of the scan.  6. Tilted 1000-particle beam in the UR model after full acceleration (left), and at beam stop (right), where the beam tilt is determined by scanning through injector phases until the one that yields the most symmetric input/output σ δ is found. Note that the center of the final distribution after the ERL is a mirror image of that before entering the accelerator. This symmetry also determines that the slope of the distribution is zero at its center after half the ERL has been traversed. Optimal beam tilt parameters for a symmetric injector-stop energy spread, found for a typical-sized CBETA beam (z = 1.1 mm, δ = 5 · 10 −4 ) with matrices (top), or by scanning through injector phases (bottom). "Angle" describes the counter-clockwise angle between the positive zaxis in longitudinal phase space and a linear fit of the beam. The beam of order 10 −3 in both coordinates is much larger than the z = 10 −7 m,δ = 10 −7 offset used to calculate the matrix, resulting in decreased accuracy for matrix tilted values.
The UR and FT models achieve the closest input and output energy spreads with injected particle distributions of similar initial tilt angles. Scan and matrix tilting meth-ods may yield different optimal linear-fit tilt angles for the same beam, because beams that have passed through an injector cavity will experience nonlinear curvature, while artificially tilted Gaussian beams will retain the uniform slope imposed in Eq. (37). For a beam with balanced acceleration and deceleration phase space profiles, injected beams that are tilted with the Table (V) angles will most effectively preserve the ERL symmetry from the optimized single-particle settings.

VII. DISCUSSION
In each of the successive models (TL, UR, FT, RK), the complexity of the RF cavity time and energy tracking has been increased successively to better model the behavior a realistic cavity. These steps were useful in determining whether the initial time symmetry conditions derived for a thin cavity situation could also be applied to finite-length, non-ultrarelativistic cavities, and finally to cavities that consider a full integrated electric field profile.
For an ideal single particle, symmetry enforcement using specific phase and flight time relations provide optimized solutions well within the tolerable range of cavity power load and peak beam energy objectives. With expected instrument fluctuations and imperfect parameter settings, solutions for the TL, UR, FT, and RK models may be more difficult to achieve, indicating a need for approximately twice the current control resolution over inputs.
Applying symmetry to a system such as CBETA allows simplification in the optimization process for ERLs like CBETA. The models in this study considered only cavity, drift pipe, and return loop pipe elements; more complex effects, such as beam optics and transverse dynamics, have not yet been examined. Further work is needed to confirm how representative the pillbox cavity models are of ERL loading and beam dynamics, as well as whether ERL symmetry is feasible for implementation in a system with elements beyond simple RF cavities and drifts. But this lies outside the scope of this paper. Nevertheless, the initial solutions found here indicate that symmetry enforcement is useful for optimizing a reduced set of phase and flight time settings to achieve cavity load and energy objectives. The phases determined here phases will be used during beam commissioning of CBETA.