Multipacting¶

Multipacting is a resonant phenomenon arising from the emission and subsequent multiplication of charged particles in accelerating radiofrequency (RF) cavities, which can limit the achievable RF power. An electron emitted from the cavity wall is accelerated by the RF field, strikes the wall again and, if the impact energy falls in the range where the surface’s secondary emission yield (SEY) exceeds unity and the field phase is synchronous, releases more than one secondary electron. Repeating over many RF periods, this avalanche absorbs RF power and can quench superconducting cavities. Predicting the field levels at which multipacting occurs is therefore crucial when optimising cavity geometries.

PyMultipact analyses multipacting in 2D axisymmetric structures in three stages, following the workflow common to multipacting codes [YlaOijala1999]:

1. solve the Maxwell eigenvalue problem (MEVP) for the resonant modes of the structure (NGSolve finite element framework),

2. integrate the relativistic Lorentz force equation to track electrons emitted from the cavity wall through the modal fields,

3. detect and resolve wall collisions, applying the SEY data to compute multipacting metrics over a sweep of field levels.

Maxwell Eigenvalue Problem¶

For the solution domain \(\Omega\), the MEVP reads

where \(\mathbf{n}\) is the surface normal, \(\mathbf{E}\) is the electric field, \(\lambda = (\omega/c)^2\) are the eigenvalues, \(\omega\) is the angular frequency, \(c\) is the speed of light and \(\partial \Omega := \partial \Omega_\mathrm{1}\cup \partial \Omega_\mathrm{2}\) is the boundary of \(\Omega\). A similar equation could be written for the \(\mathbf{H}\) field. The domain \(\Omega\) is assumed to be vacuum.

The unknown field is expanded in basis functions of a finite element space,

where the natural functional space for the \(\operatorname{curl}\operatorname{curl}\) problem is \(H(\operatorname{curl}, \Omega)\); the \(\mathbf{u}_i\) are the edge (Nédélec) basis functions.

For an axisymmetric domain, the volume integrals of the weak (Galerkin) form pick up the radial coordinate \(r\) from the cylindrical volume element \(\mathrm{d}V = r\,\mathrm{d}r\,\mathrm{d}z\,\mathrm{d}\varphi\), so the weak form solved by PyMultipact on the 2D \((z, r)\) cross-section is the r-weighted problem

where \(\mathbf{v}\) is the test function (equal to the basis functions). In matrix form this is the generalised eigenvalue problem

PyMultipact assembles both forms with NGSolve and solves the eigenproblem with the preconditioned inverse iteration (PINVIT) solver. Because the \(\operatorname{curl}\) operator has a large kernel (all gradient fields), a divergence-free projector built from the discrete gradient matrix is applied inside the preconditioner to suppress the spurious zero-frequency space and keep the iteration on the physical modes.

Third-order Nédélec elements are used by default (Domain.compute_fields(order=3)), matching the element order of the MultiPac field solver [YlaOijala2001]; this retains field accuracy on a comparatively coarse mesh. The boundary polyline of dense geometry files is resampled to ~250 points by default (Domain(n_boundary_points=...)) so the mesh size follows the element order rather than the input file resolution.

Relativistic Lorentz Force¶

The motion of charged particles in electromagnetic fields is described by the Lorentz force equation. Since the emitted electrons are relativistic (\(\beta = 1\)) or near relativistic (\(\beta \approx 1\)), the relativistic form [YlaOijala1999] is integrated:

where \(q\) and \(m\) are the charge and mass of the electron.

The launch ensemble is the phase space of surface position and RF phase,

i.e. each initial particle is a pair \((\mathbf{x}_i, \psi_i)\) of an emission site on the wall and an initial field phase. Particles are released with a small initial velocity \(v_0\,\mathbf{n}\) along the inward surface normal (\(v_0\) corresponding to e.g. 2 eV) and tracked for a fixed number of RF periods at every field level of the sweep.

Integration scheme¶

Classic Runge-Kutta¶

With the right-hand side

the classical fourth-order Runge-Kutta scheme is used. Let the values at time \(t_n\) be \(\mathbf{v}_n\), \(\mathbf{x}_n\) and the time step \(h\):

If an intermediate stage carries a particle outside the domain, the stage is resolved by the collision handler (below) before the update is completed.

Collision Detection and Resolution¶

Let

be the line segment connecting the particle positions at time steps \(k\) and \(k+1\), and write the discretised boundary as \(\partial\Omega := \bigcup_i \partial_i\Omega\) of straight wall segments. A collision occurs iff \(|L_k \cap \partial_i\Omega| > 0\) and odd. The candidate wall segments are found with a k-d tree over the boundary points; the intersection point \(\mathbf{x}_\mathrm{c}\) and the exact impact time within the step are recovered from the crossing.

At the moment of impact the electric field at \(\mathbf{x}_\mathrm{c}\) is evaluated (from the FEM solution). If the field pulls electrons away from the wall, i.e. \(\mathbf{n}\cdot\mathbf{E} > 0\) with \(\mathbf{n}\) the inward normal, a secondary is released with the same initial speed \(v_0\,\mathbf{n}\) and tracking continues; otherwise the particle is counted as absorbed (loss_model='field', the default and the behaviour benchmarked in the paper). Two alternative loss models can be selected in Domain.analyse_multipacting(loss_model=...) for numerical experiments: 'wait' re-emits the electron uncounted until the RF phase turns favourable (similar in spirit to delayed re-emission), and 'always' re-emits and counts every impact.

Multipacting Metrics¶

The counter function (CF) and enhanced counter function (ECF) of the Helsinki group [YlaOijala1999] are implemented.

Counter function. \(c_n/c_0\) is the ratio of the number of particles that survive at least \(n\) impacts (by default \(n=20\)) to the number of particles released. It identifies surface sites and field levels that geometrically support multipacting.

Impact energy. The kinetic energy at the \(k\)-th impact is

with \(\gamma_k\) the Lorentz factor at impact and \(q_0\) the elementary charge. The final impact energy \(Ef_{20}\) reported by Domain.plot_Ef is the mean of the last (20th) impact energies of the surviving particles; multipacting is dangerous where \(Ef_{20}\) falls in the region where the SEY exceeds unity (between the crossover energies of the SEY curve, drawn as reference lines).

Enhanced counter function. The number of secondaries generated by one particle after \(n\) impacts is

where \(\delta\) is the SEY of the surface (interpolated linearly from tabulated data, e.g. measurements). For \(P\) initial particles the ECF is

\(e_{20}/c_0 > 1\) indicates that the electron population actually grows, i.e. genuinely dangerous multipacting rather than merely resonant trajectories.

Distance function. For each surviving trajectory the distance in (position, phase) space between the final and initial state,

measures how exactly the electron returns to its starting point: small \(d_{20}\) identifies stable (fixed-point) multipacting orbits. Domain.plot_df renders it as a map over emission site and initial phase.

For two-point multipacting the orbit returns to the launch side only on every other impact, so the plotted default compares the launch point with the nearest of the last two impacts — this removes the arbitrary branch parity of the impact count and reproduces MultiPac’s graded map (dark core at the resonant fixed-point phase, growing with the launch offset). The literal 20th-impact definition (metric='d20_strict') and the closure of the two-impact map (metric='closure', zero for any phase-locked orbit) are also available.

References¶