{'fraction_non_alphanumeric': 0.040721912434695424, 'fraction_numerical': 0.015532375239670004, 'mean_word_length': 4.371315192743764, 'pattern_counts': {'":': 0, '<': 5, '<?xml version=': 0, '>': 5, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 1, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 10, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Field Analysis for a Highly-Overmoded Iris Line, with Application to THz Radiation Transport at LCLS

Vector field analysis is presented for a highly overmoded iris-line structure that can act as a medium of transportation for THz radiation. The axisymmetric structure is capable of supporting hybrid modes with desirable features such as low propagation loss, uniformly linear polarization and approximately-Gaussian intensity profile across the iris. A specific application that can benefit from these features is the transportation of THz undulator radiation over hundreds of meters to reach the experimental halls at the LCLS facility, SLAC, Stanford. Such a structure has been modelled before as a boundary-value problem using Vainstein's impedance boundary condition and assuming infinitely-thin screens. Given that physical realizations of such screens must have finite thickness and that the THz wavelength in the 3-15 THz range is expected to be smaller than convenient screen thicknesses in practice, the question of the impact of finite screen thickness on propagation performance becomes rather pressing. To address this question, we present a mode-matching analysis of the structure, with finite screen thickness and perturbatively clustered (localized) field expansions, for computational feasibility. The effect of screen thickness is seen to lower the attenuation constant on the iris line, which is dominated by diffraction loss. The propagation loss predictions based on the Vainstein-model are compared with the numerical results from mode-matching for infinitely-thin screens, where the former method is observed to agree with the numerical results better at higher Fresnel numbers (highly-overmoded structures). The properties of the dominant mode fields are formally derived from first principles and a recommended approach is discussed for the inclusion of screen-thickness effects into propagation loss estimations.

Vector field analysis is presented for a highly overmoded iris-line structure that can act as a medium of transportation for THz radiation. The axisymmetric structure is capable of supporting hybrid modes with desirable features such as low propagation loss, uniformly linear polarization and approximately-Gaussian intensity profile across the iris. A specific application that can benefit from these desireable features is the transportation of THz undulator radiation over hundreds of meters to reach the experimental halls at the LCLS facility at SLAC, Stanford. Such a structure has been modelled before as a boundary-value problem using Vainstein's complex-impedance boundary condition and assuming infinitely-thin screens. Given that physical realizations of such screens must have finite thickness and that the THz wavelength in the 3-15 THz range is expected to be smaller than convenient screen thicknesses in practice, the question of the impact of finite screen thickness on propagation performance becomes rather pressing. To address this question, we present a mode-matching analysis of the structure as an open resonator with finite screen thickness and perturbatively clustered (localized) field expansions, for computational feasibility. The effect of screen thickness is seen to lower the attenuation constant on the iris line, which is dominated by diffraction loss. Ohmic loss due to the finite conductivity of metallic surfaces at the screen edges are found to be negligible compared to diffraction loss. The propagation loss predictions based on the Vainstein-model are compared with the numerical results from mode-matching for infinitely-thin screens, where the former method is observed to agree with the numerical results better at higher Fresnel numbers (highly-overmoded structures). The properties of the dominant mode fields are formally derived from first principles and a recommended approach is discussed for the inclusion of screen-thickness effects into propagation loss estimations.

I. INTRODUCTION
Efficient transportation of THz radiation over long distances (hundreds of meters) is a challenging problem. On one hand, the THz wavelength (say, 20-100 µm, for the range 3-15 THz), being smaller than microwaves, makes the scaled transverse dimensions of a traditional singlemoded waveguide impractically small over long distances. On the other hand, being larger than optical wavelengths makes the THz wave relatively more prone to Fresnel diffraction. An example of this challenge, which was the original motivator behind the present study, is the problem of efficiently transporting radiation from an "afterburner" THz linear undulator downstream of LCLS over a distance of 150-350 m, to reach the experimental halls at the LCLS facility, SLAC, Stanford [1,2]. A traditional quasi-optical solution that utilizes a combination of planar, toroidal or paraboloidal mirrors to relay the THz beam in multiple steps is one proposed solution, which typically suffers from power loss of approximately 1% per mirror as well as some aberration and misalignment [1,3]. To reach the near experimental hall at LCLS, for example, through a 150-m path (roughly 34 mirrors) going through the access maze at LCLS, the mirrors are estimated to incur around 30% power loss [1]. The iris-line structure, whose analysis is the subject of this paper, is an alternative solution that was first proposed by Geloni, * anaji@stanford.edu; et al., [2] for the THz transport at LCLS. The iris line supports an attractive hybrid mode that can deliver (1) a low propagation loss, (2) an almost-Gaussian intensity profile and (3) an invariant linear polarization across the iris. These features make it ideal for direct coupling with the radiation incident from the THz linear undulator.
The iris-line geometry under consideration consists of a series of parallel screens, each having a circular iris (gap) of radius a and centered around the axis z, as shown in Figure 1. The periodic structure has a period b and iris radius a that are much larger (hundreds or throusands) of the THz wavelength in vacuum, λ 0 , making the structure highly overmoded. The screens expand transversely to an outer radius r 0 > a that defines the boundary of the enclosing chamber. The "virtual pipe" furnished by the irises along the axial (z) direction provides the resonant structure with a means to act as a waveguide (transmission line). On the other hand, the structure's resonance conditions, controlled by the periodicity and the shape of the annular regions sandwiched between the screens, will constitute equivalent boundary conditions from the perspective of the waveguide.
For a plane wave that is paraxially incident on the iris line (Figure 1), the Fresnel diffraction (or knife-edge diffraction) experienced at the screen edges will cause part of the wave to be lost into the shadow region between screens. For a single screen edge illuminated by a plane wave, this is similar to the well-known phenomenon of knife-edge obstruction (e.g. [4]), which causes the diffracted field amplitude to decrease as we go deeper into the shadow (see Figure 2a) according to the atten- uation factor |F (ν)|, where F (ν) = 1+i 2 ∞ ν e −iπτ 2 /2 dτ and ν is the Fresnel diffraction parameter defined as ν = h 2(d 1 + d 2 )/(λ 0 d 1 d 2 ). Assuming that the screens are perfectly absorptive, all the diffracted waves entering the shadow will constitute lost power from the perspective of propagation along the axis. This can be equivalently modeled as an open-resonator structure, where we let the outer radius r 0 go to infinity and the screens can be assumed to be perfectly absorptive or conductive; in either case the same diffraction mechanism is seen at the screen edges and none of the diffracted waves entering the shadow region will return to the line, as shown Figures 2b,c. Note that such an ideal open-resonator model (i.e. with r 0 →∞) can be well approximated in practice by making the distance r 0 −a finite but deep enough such that the power lost to diffraction is approximately equal to that lost in the limit r 0 →∞. This is typically achieved by taking the depth r 0 −a to be much larger than the diffraction scale ∼ √ bλ 0 , [5]. One way to analyze this open resonator is to consider the diffraction effects between the screens equivalent to a complex impedance boundary condition at r = a and  2. (a) A canonical example of Fresnel diffraction of a plane wave caused by a single knife-edge obstruction (thin and absorptive screen). The field strength of the wave traveling from point S to point P diminish according to the function |F (ν)| shown. As the point P moves deeper into the shadow region behind the screen, the attenuation is increased. The Fresnel parameter ν = h 2(d1 + d2)/(λ0d1d2) defining the depth of shadow (or clearance) is found geometrically from the dimensions shown. (b) A simplified illustration of rays lost to diffraction, as they hit an absorptive screen or get diffracted into the shadow region bound by two such screens. (c) A simplified illustration of rays lost to reflection or diffraction in the region between two conductive screens. thereby convert the problem into a boundary-value problem over the closed virtual pipe of radius r = a, [2]. Such a boundary condition has been called the Vainstein boundary condition [2,6,7] and provides one of the methods of analysis we shall use to describe the losses on the iris line. This method, however, assumes infinitely-thin screens (δ = 0). In practice, the screens will have to be of finite thickness. Indeed, the finite thickness of the screen will be quite large relative to the wavelength at THz frequencies and a question is raised as to whether the screen thickness will influence the attractive features (low loss, invariant linear polarization and almost-Gaussian amplitude profile) of propagation using the iris line, as predicted by Vainstein's model. To address this question we analyze the open-resonator structure using the mode matching method and clustered (localized) field expansions, for computational feasibility.
For this type of overmoded open-resonator structures, it will be seen that the dominant power loss mechanism is due to diffraction or radiation loss, not ohmic (conductive) loss. This should be compared with smooth waveguides and traditional corrugated structures, typically found at microwave frequencies where the depth and period of corrugations are shallow and often ≤ λ/4 (e.g. [8][9][10][11][12][13]). For example, the dominant hybrid mode on the iris line will be shown to have the desirable property of the attenuation constant decreasing with increasing radius and increasing frequency, as a function of a −3 ω −3/2 . The same property for the attenuation constant is found in smooth circular waveguides operating with the TE 0n modes [8]. However, the attenuation on the iris line is mainly due to diffraction loss, while on the smooth circular waveguide is due to ohmic loss (by azimuthal surface currents) [8,14].
Employing this iris line for THz radiation transport in practice will clearly require a study of not only its propagation power loss, polarization purity and amplitude profile, but also of critical considerations such as correct mode launching and coupling, input transients, dispersion and tolerance to fabrication errors and mechanical misalignment. In this paper, we confine our scope to the investigation of propagation loss, polarization and amplitude profile properties, assuming screens of finite thickness and paraxial plane wave incidence (no electron beam on the line).
The paper is organized as follows. In Section II we briefly review the Vainstein-based model and use its predictions as a reference for comparison. We then proceed to Section III where we derive the vector field equations for the dominant dipole mode of the iris line using the method of mode matching and clustered expansions. In Section IV we discuss the implementation of modematching and present its results for the specific example of the iris-line structure proposed for THz transport at LCLS, while investgiating the effect of finite screen thickness on the propagation properties. The paper concludes in Section V, followed by two Appendices that contain more formal derivations. Appendix A presents the use of perturbation theory to formally derive the vector field properties of the hybrid mode on the iris line starting from Vainstein's model. Appendix B uses perturbation theory to provide a justification for the use of the method of clustered field expansions in the analysis of overmoded paraxial iris lines.
Throughout this paper, we assume harmonic time dependence of the form e −iωt .

II. ANALYSIS BASED ON VAINSTEIN'S COMPLEX BOUNDARY CONDITION
The open-resonator model of the iris line can be converted to an equivalent closed resonator problem which consists of Vainstein's approximate boundary condition at the cylindrical wall r = a. This boundary condition is formulated by Geloni,et al.,in [2] as being the Fresnel number, and E represents the electric field envelope components E r or E θ in the cylindrical frame of coordinates (r, θ, z). This analytical model assumes that the parameter M is small (M 1) and that the structure is highly-oversized compared to the wavelength, with k 0 b 1, where k 0 = 2π/λ 0 is the wavenumber in free-space.
Solving the Helmholtz wave equation while imposing the boundary condition (1) results in a dominant hybrid mode with invariant linear polarization across the iris, an amplitude profile approximately equal to the Bessel function J 0 (2.4r/a), and the following propagation power loss (L p ) law (see Appendix A for derivations) where c is the speed of light in vacuum, ω is the angular frequency and z is the distance travelled down the line. A key feature of this loss law is the tendency of the attenuation constant Im[β] to drop with larger radii and frequencies (inversely proportional to a 3 and ω 3/2 ), where β is the complex propagation constant along the z direction. An iris line with dimensions a = 5.5 cm and b = 30 cm, for example, has been proposed [2] for the THz radiation transport at LCLS. Using the power loss estimation (2), we see that a THz undulator radiation in the range 3-15 THz will experience maximum propagation loss of 14% at 3 THz. The loss decreases at higher frequencies, as shown in Figure 3. If larger radii are feasible for installation at LCLS, lower loss can be achieved by exploiting the cubic-law dependence of loss on radius a. For example, Figure 3 shows a second iris line with a = 10 cm, and a power loss less than 3% at 3 THz.
The predictions given by (2) do not take into account the effect of finite screen thickness (δ > 0). In Section III, we analyse the modal fields on the line for nonzero screen thicknesses and compare the results with those predicted by Vainstein's model.

III. ANALYSIS BASED ON MODE MATCHING IN AN OPEN RESONATOR
In this analysis we begin by expanding the modal fields in the axial region (region 1 in Figure 4) and the annular region (region 2 in Figure 4), then impose the boundary conditions at r = a and derive the structure's characteristic dispersion relation (k 0 -β) in matrix form. Solving the characteristic equation for a given wavenumber k 0 = ω/c will then yield the complex propagation constant β on the iris-line. The presented approach is similar to the classical treatment given by Zotter and Bane [15] in terms of matrix equation conditioning. However, the present treatment is different from [15] in that it is concerned with the analysis of an open (r 0 →∞), rather than closed, resonator model and in that we derive all the field components directly using the longitudinal components (E z , H z ) of the fields themselves, rather than the Hertz potentials. Recall from Whittaker's theorem [16,17], that only two independent solutions of the scalar wave equation are needed to determine all the 6 components of the EM field in vacuum. The two approaches are therefore fundamentally equivalent.
A. Field expansion in the axial region (region 1) Starting with Maxwell's curl equations, ∇×E = iωµH and ∇×H = −iω E, and seeking wave solutions that are travelling along z in the form e iβz , we can rewrite these curl equations with all the vectors and the operator ∇ divided into their axial (along z) and transversal (subscripted with t) parts, as For a transverse electric (TE) mode, we substitute E z = 0 into (3) and (4) and equate vector terms that are parallel to each other on each side of the equations. After algebraic vector manipulation, we obtain For a transverse magnetic (TM) mode (H z = 0), we similarly find If the boundary conditions of a given structure provide a form of coupling (impedance) that links the TE and TM modes, a hybrid mode is generally produced as a combination of the TE and TM modes, with Consider now our iris-line structure and the field description for region 1, which must generally meet the hybrid mode equations (9) and (10), as well as the scalar transverse Helmholtz equations, ∇ 2 t E z + k 2 t E z = 0 and ∇ 2 t H z + k 2 t H z = 0, where k t = k 2 0 − β 2 is the transverse wavenumber. Since both the Helmholtz equation and the present boundary conditions are separable in cylindrical coordinates [18], we can immediately assume that E z and H z are of the form E z = CR(r)Θ(θ)Z(z) and H z = DR(r)Θ(θ)Ẑ(z), where C and D are arbitrary constants. We find the form of the functions R(r), Θ(θ) and Z(z) that constitute E z as follows (with a similar treatment for H z ). Given the periodicity of the structure, Z(z) is predicted by Floquet's theorem to have the form of a Bloch wave [19,20], whose field will be periodic from one structural period to the next, except for a complex phase advance of e iβ0b , where β 0 is the propagation constant along z. Note that β 0 is generally complex to allow for attenuation in lossy structures (as in our iris-line structure), but would be real when the structure exhibit no losses to materials or radiation (i.e. closed structure with perfectly conducting walls and perfect dielectrics) [20]. We can therefore write Z(z) as Z(z) = e iβ0z φ b (z), where φ b (z) is a periodic function of fundamental period b. Expanding φ(z) as a spacial Fourier series, n C n e i2nπz/b , leads to where β n is the propagation constant and k tn = k 2 0 − β 2 n is the transverse wavenumber for the n th spacial harmonic. As in (9) and (10), the transverse components of the n th harmonic can be found from its E zn , H zn components as To find the radial and azimuthal dependences, R and Θ, we now write the Helmholtz equation ∂θ 2 ] in such coordinates couples the radial and azimuthal components without coupling the axial component, the Z(z) function drops out of the equation and we end up with the following equation, after separating the variables, where m is a constant integer (to keep the azimuthal dependence a single-valued function of θ). This gives azimuthal dependence Θ in the form of cos mθ or sin mθ. With no loss of generality, we use the former for E zn and the latter for H zn , since the azimuthal dependenceΘ in H z will turn out (as it must) to be of the same form.
We also anticipate that m will be the same across region 1 and 2, to maintain phase matching and as it does not depend on n; we hence reuse the same symbol m in both regions. The radial dependence is now seen to reduce to which is the parameterized Bessel equation. Since the radius in region 1 is bounded by the iris radius (r ≤ a) and includes the axis (r = 0), we choose our solution to be the Bessel function of the first kind, J m (k tn r), and reject the second kind. With a similar treatment for H zn , we now have where, for symmetry, we have chosen to give C and D the same physical dimensions and explicitly isolate the free-space impedance, Z 0 = µ/ . We now substitute (16) and (17) into (12) and (13) to yield, after algebraic manipulation, the hybrid field expressions for region 1, as where J m [·] denotes the derivative of J m [·] with respect to its argument. Note that these equations will be further normalized following mode-matching in Section III C, to provide symmetry and convenience during calculation.

B. Field expansion in the annular region (region 2)
Since the field propagation in region 1 will vary periodically along z with the same fundamental period as the structure's, we know that the field in one annular gap will be the same as that in all gaps. Therefore, we merely need to solve for one gap. Since the gap is bounded in the z direction and open in the r direction, we expect it to host a standing wave along the longitudinal direction and a traveling wave in the radial direction (radiating outwards). The latter is expected to meet the usual radiation boundary condition (E, H → 0) at infinity. For the fields E z , H z in the gap, we also expect to superimpose solutions of the form e ±iβpz in order to form the standing waves. These fields must satisfy the transverse Helmholtz wave equations, p and the index p denotes a given mode (standing wave) in region 2. Since both the Helmholtz equation and the structure's boundary conditions in region 2 are separable in the cylindrical frame [18], we follow a treatment similar to the one used above for region 1 and we assume that For the p th TE mode, considering that H z must have nodes at the screen walls (see Figure 4), thenẐ p (z) must be of the formẐ Similarly, for the p th TM mode, considering that E z must have maxima at the screen walls, Z p (z) must be of the form . We now note that equations (9) and (10), which related the transverse components to the longitudinal fields (E z , H z ), were originally derived for a travelling wave (waveguide setup) rather than a standing wave (gap setup), and will need to be adjusted accordingly. For the TE (E z = 0) case, we now need to take ∂/∂z ≡ iβ p for the e iβp(z+∆) term and ∂/∂z ≡ −iβ p for the e −iβp(z+∆) term; and similarly for the TM mode. Repeating the derivation steps taken for region 1 by separating Maxwell's curl equations in terms of the longitudinal and transverse parts for the p th mode, followed by algebraic vector manipulation, gives the following relations for the TE mode In a similar manner, we have for the TM mode To find the radial and azimuthal functions (R, Θ) in this region, we substitute into the Helmholtz equation and use the cylindrical Laplacian, as was done in region 1. This leads to the same azimuthal dependence (namely, cos mθ or sin mθ) and the same parameterized Bessel equation (15). For its Bessel solutions in region 2, we pick the Hankel function of the first kind, H m at large distances [21], viz.
Substituting the R, Θ,R,Θ forms back into equations (24)-(29), we can now explicitly write the full hybrid field expressions for region 2 as where B 0 = 0. Note that, for symmetry in equations and parameters, we take the impedance Z 0 explicitly out of the B p coefficients. These equations will be further normalized following mode-matching in Section III C, to provide symmetry and convenience during calculation.
It is important to bear in mind that the fields (31)-(36) are givens for |z| ≤ ∆ inside the gap and are identically 0 in the region ∆ < |z| ≤ ∆ + δ/2, where the screens are assumed to be made of a perfect conductor (see Figure 4).

C. Mode matching at r = a
We now solve for the coefficient families (A p , B p ) and (C n , D n ) by enforcing the continuity of the tangential E and H fields across the boundary r = a (the continuity of the normal fields follows automatically from Maxwell's divergence equations). Namely, we require the four conditions at r = a Enforcing the first condition gives us Multiplying both sides by e −iβnz and integrating through 1 b b/2 −b/2 , using integration by parts twice, this yields (after manipulation and abbreviating the obvious Bessel argument) equations (39) and (40) below. Following similar steps, we enforce the second, third and fourth conditions, using integration by parts twice for each case, to also get equations (41)-(43).
sin β n ∆, for even p, −i cos β n ∆, for odd p, sin β n ∆, for even p, −i cos β n ∆, for odd p, sin β n ∆, for even p, −i cos β n ∆, for odd p, sin β n ∆, for even p, i cos β n ∆, for odd p, where the dagger symbol represents Hermitian conjugation and δ p0 is Kronecker's delta (equal to 1 when p = 0 and zero otherwise). Note that M and N will denote matrices that contain elements such as M np or N np .
Further simplification and symmetry in the equations are achieved if we momentarily define B p ≡ iB p and remove the imaginary i factor from odd cases in the M, N matrices by absorbing it into the correspondingĀ p ,B p coefficients, calling the new matrixedM ,N , as follows sin β n ∆, for even p cos β n ∆, for odd p (46) sin β n ∆, for even p cos β n ∆, for odd p This would render the equations real, since the coefficients themselves are arbitrary. For symmetry and convenience, we also normalize the coefficients as which is equivalent to dividing the original mode equations in (18) and (31) Using matrix conditioning and reduction steps similar to those used in [15], we now take the following definitions to simplify our equations, where all the matrices in (60)-(64) are diagonal (either in p or in n) and, therefore, easy to manipulate. The modal equations can now be written concisely as In full matrix notation we can cast this into whereĀ,B, C, D are vector matrices,N ,M are full real matrices and the remaining matrices are diagonal. We can further write this system of equations in terms of "matrices of matrices" as follows where 1 and 0 here represent the unity and zero matrices. We can solve these equations to eliminate theĀ,B pair. After manipulation, we obtain By appropriate grouping of the matrices in this homogeneous result we can finally write the system's characteristic equation explicitly in symmetric form as whose determinant must vanish to give the eigensolutions and the dispersion relation between ω and β 0 , after we truncate the p sum to some index P max and the n sum to ±N max in practice (remembering that the index n runs in region 1, whilst the index p runs in region 2). Specifically, we require that sin β n ∆, for even p −i cos β n ∆, for odd p (83) sin β n ∆, for even p −i cos β n ∆, for odd p Solving (77) numerically for a given frequency or wavenumber k 0 = ω/c will produce the sought complex propagation constant β 0 of the iris-line, whose real part Re[β 0 ] will represent the longitudinal phase constant and imaginary part Im[β 0 ] the longitudinal attenuation constant. Once we have found β 0 , we can write the full fields by substituting into (48)-(59).

IV. NUMERICAL IMPLEMENTATION AND RESULTS
The mode-matching model obtained above is now ready for numerical implementation on a computer. Our investigation of the iris line is concerned not only in comparing the numerical results with the predictions of the model based on Vainstein's boundary condition, but also observing the influence of finite-screen thickness on the attractive characteristics (propagation loss, polarization and amplitude profile) of the dominant mode. For the open iris-line structure in hand, we expect radiation loss and the paraxial propagation (similar to a plane wave, but with slightly bent, paraboloidal wavefronts [22]) to contribute, respectively, by giving a nonzero imaginary part of β 0 and a real part that is close but not exactly equal to k 0 [see equation (B9) in Appendix B for a more formal description].
Assuming an incident THz wave with linear polarization, we are mainly interested in the dominant dipole mode on the iris line that can couple to and support such a wave. The basic aim of the numerical code is therefore to find the zeros of the determinant (77) for a field with m = 1, which would give us β 0 = Re[β 0 ] + iIm[β 0 ]. Given the periodicity of the dispersion curves in such a system, there will be an infinite number of possible solutions. The dominant mode will naturally be the one with the lowest attenuation constant (smallest Im[β 0 ]). Given the homogeneous equation (76) of the C n , D n field coefficients, our requirement for the determinant of this matrix to vanish is basically a requirement to find a value of β 0 that gives nontrivial field solutions. Substituting such a β 0 back in to the matrix will thus lower its rank and link the field coefficients C n , D n to each other, informing us about the "shape" (or profile) of the dominant vector field distribution. Solving this in practice using standard computer tools, such as Wolfram Mathematica (v.12.3.0), can be formulated as an eigenvalue problem. Indeed, the characteristic matrix (76) and [C, D] T vector can be written in shorthand notation as Ax = 0, where A denotes the matrix (76), x the [C, D] T vector, and the superscript T denotes taking the transpose. Then solving the homogeneous equation Ax = 0 is equivalent to solving the eigenequation Au = λ u u to find the eigenvector u that corresponds to the zero eigenvalue (i.e. if λ u = 0, then x ≡ u). Taking the C n and D n coefficients from this eigenvector and subtituting them, alongside β 0 , into (48)-(59) will then allow us to visualize the vector field polarization and profile, as will be shown in Subsection IV C.
In the following subsections, we discuss the various aspects related to physical interpretation and practical implementation of the mode-matching method.
A. Discussion on the truncation and clustering of field expansions For numerical implementation, one clearly needs to truncate the infinite series in the n and p expansions seen in (48)-(59). Ideally, the expansion in n would be , where the value of P max is customarily chosen close to 2N max (to facilitate matching both expansions across boundary, assuming thin screens) and both are high enough to allow for the solution to numerically settle at a given accuracy. Indeed, the value of P max can be estimated by intuitively considering how the field in region 2 is excited by the wave incident in region 1. Since we assume paraxial incidence (β 0 ≈ k 0 ) from the THz source along the z axis, and since the proposed structure is highly overmoded, with each period b in the order of ∼3000 wavelengths, the longitudinal (guided) wavelength will roughly be equal to its value in free-space ∼λ 0 = ω/c. When this wave is matching tangentially across the boundary (r = a) to- In this sketch βn is assumed to be real (no loss). For a smooth pipe, the continuous parabola (green) is observed, where phase velocity v ≥ c and k0 ≥ β0 for propagating modes, asymptotically leading to the limit of TEM propagation at the 45-deg line. As we start loading the pipe with shallow (weak) irises, the propagated modes are allowed in certain passbands, attenuated in stopbands, and the dispersion curves (red) become periodic in βn, with every side lobe period (index n) now corresponding to a Fourier harmonic bin, where βn = β0 + 2πn/b. This allows for slow and fast waves to propagate (depending on the used harmonic). If the irises are gradually removed, the passband curves (red) expand and coalesce back into the parabolic curves of the smooth pipe (green) [20]. If the iris loading increases, the passbands (blue) are compressed, and in the limit they become flat spectral lines, with wider stopbands in-between. Since our structure is highly overmoded with a paraxial incidence, we effectively operate in the perturbation region (purple shade) in the vicinity of the parabolic curve and make our choice for the n = 0 harmonic band accordingly; see the simplified example of point (β * 0 , k * 0 ) on the figure. wards region 2, whose width is 2∆ (see Figure 4), we expect the excitation of the mode whose longitudinal standing-wave in region 2 has a half-wavelength value close to the value of ∼λ 0 /2. This implies that the dominant mode number (call it p = P 0 ) in region 2 will be approximately the integer closest to 2∆ λ0/2 . This, indeed, turns out to be the case, as can be observed through numerical computations.
The modes in region 2, however, can be categorized, as in traditional smooth waveguides, to either propagating or evanescent modes. The former happens when k 0 > β p (giving real k tp , with phase velocity v > c and H (1) m (k tp r) Hankel function radiation), while the latter corresponds to k 0 < β p (giving purely imaginary k tp = iγ p , where γ p is a real number, with v < c and K m (γ p r) modified Bessel function radiation). The dispersion curve of such a waveguide takes the usual form of a parabola that is entirely above the 45-deg line (speed of light line), asymptotically approaching it from above in the high-frequency limit; see Figure 5. Therefore, the integer P 0 closest to 2∆ λ0/2 must actually be chosen as the floor integer (i.e. the closest integer from below) to represent a propagating mode (radiation) in region 2. Such a mode will typically be responsible for most of the attenuation seen in Im[β 0 ]. If the ceiling integer was chosen instead, it would correspond to the largest evanescent mode. Therefore, in our implementation and throughout the results below P 0 will refer to the floor integer of 2∆ λ0/2 = 4∆ λ0 . The choice of P 0 also serves as a natural midpoint around which we may choose to "cluster" (or localize) our p expansion, so that, instead of summing over Pmax p=0 , we sum over P0+pSteps P0−pSteps , where we choose pSteps during numerical implementation. A tradeoff clearly exists between higher accuracy of the results (higher pSteps) and the computational burden (time and memory) required to perform the calculation. Our choice to use, or not use, clustering will depend on how overmoded a structure is and how difficult it is to run a regular expansion (without clustering). Note that a regular expansion Pmax 0 is equivalent to a clustered expansion P0+pSteps P0−pSteps , if pSteps happened to be equal to P 0 and P max happened to be equal to 2P 0 . Moreover, if a regular expansion up to P max stopped short of reaching P 0 (that is, if P max < P 0 ), no correct solution will be found, because the dominant mode would not be included. This highlights the usefulness of clustered expansion when we can only expand in a limited number of terms (e.g. due to computational limitations), since the clustered expansion will always include the dominant mode; see Figure 6.
A similar concept can be used to cluster the n expansion around the dominant harmonic, which is the closest (rounded) integer (call it n = N 0 ) to the expected value b λ0 . Unlike P 0 , there is no need to seek the floor or ceiling integer values for N 0 . Furthermore, we also need to include the Fourier image of the dominant harmonic.
At this point, we are faced with a question as to where to look for the value of β 0 on the dispersion curves, when such curves are periodic in β n = β 0 + 2πn/b? (see Figure 5). The periodicity naturally implies that all harmonic bins (lobes) are equivalent and we can choose to designate the index n = 0 to any of them. However, since our structure is a highly-overmoded structure with a paraxial wavefront that can be considered a small perturbation compared to the TEM limit on the dispersion diagram, the dispersion curves forming the passbands in Figure 5 will be in the vicinity of the parabolic curve of a smooth waveguide (small iris-loading effect). Indeed, in the limit of no iris loading (no diffraction losses), these passbands coalesce back to give the parabolic curve [20]. For a perturbative iris effect, the dominant harmonic is expected to be near the corresponding (k 0 , β 0 ) point on the parabolic curve in Figure 5. Consequently, we choose n = 0 to be the index for the dominant mode corresponding to β n = β 0 ≈ k 0 , which is equivalent to N 0 ≈ b/λ 0 harmonic bins away from zero on the β n axis. The Fourier image of the dominant mode will hence be at n = −2N 0 , and a clustering of the form N0+nSteps N0−nSteps + −2N0+nSteps −2N0−nSteps may be used for practical computations.
The advantage of working around a dominant mode that is slightly perturbed (i.e. clustering the field expansion) seems to, in some sense, offset the computational disadvantage of operating with such a highly-overmoded structure. Had the iris loading been too strong [e.g. small Fresnel number N f = a 2 /(bλ)], more wave reflections between the irises would have resulted in several strong harmonic terms, which are not necessarily the original one near the parabolic curve (the periodic dispersion curves in Figure 5 would be relatively flattened) and a preferred designation for the n = 0 harmonic bin would no longer be clear. Good discussions on the interpretation of the n harmonic bins and the limiting cases can be found in references [20] and [11].
The discussion above introduces the method of clustering in light of its intuitive nature in practice and the behavior infered from the k 0 − β dispersion curves. A more formal justification for using clustering in highly-TABLE I. The three different scales considered during the present analysis: mildly-overmoded (Scale-1), overmoded (Scale-2) and highly-overmoded (Scale-3). Scale-3 is the structure proposed for the LCLS THz transport problem. For each scale, N0 indicates the structure period (b) measured in wavelengths, N f =a 2 /(bλ0) is the Fresnel number and the values given for screen-thickness δ are example values to indicate the scaling. All shown dimensions are in mm. overmoded paraxial lines is provided in Appendix B, for completeness.

B. Structure scales
Since the iris-line dimensions (a, b) under consideration are in the order of hundreds or thousands of wavelengths at THz frequencies, making the line highly overmoded, the computational cost experienced for mode-matching analysis and iterative optimization can be prohibitively high. To help develop our understanding of overmoded iris-line structures and to help optimize our numerical algorithms, we utilize the clustering method where possible and we solve the problem over 3 scales. While keeping the aspect ratios fixed, each scale represents a similar structure with dimensions scaled 10-times relative to those of the next scale. Scale-1 dimensions are 100-times smaller than our intended structure (thus, 100 times less overmoded), Scale-2 is 10-times smaller than our intended structure, and Scale-3 is our intended structure; see Table I. All scales are analyzed at the same frequency of 3 THz. Clearly, solving the Scales 1 and 2 are less computationally prohibitive than Scale-3, where the structure period span thousands of wavelengths. Working over scales will also allow us to spot different behavior trends as the scale moves from mildly-overmoded to highly-overmoded, as will be reported in Subsection IV C. When validating Vainstein's model using mode-matching, it is important to keep an eye on the small parameter M and the number k 0 b in each scale, since Vainstein's model assumes M 1 and k 0 b 1, which imply a better validity in highly-overmoded structures. Thus, Scale-3 is expected to have stronger relevance to the predictions by Vainstein's model, compared to Scale-1 and Scale-2.

C. Numerical results
The results of mode-matching calculations are summarized in Table II for the three scales of the iris line listed in Table I. In addition to comparing the mode-   Table II for details.  Figure 8 that the full-scale structure (Scale-3) will have 33.9% reduction in the attenuation constant Im[β 0 ] for a screen that is 2-mm thick. The change is less acute for Scale-2 structures around equivalent relative-thickness points (δ/b). Scale-1 exhibits some oscillatory behavior for Im[β 0 ] as a function of δ for relatively thin screens, before it settles in value for thicker screens. This may be attributed to the small Fresnel number of Scale-1, effectively taking the structure from the perturbative diffraction limit into stronger iris loading (admittance) where more harmonic terms experience stronger reflections and can interfere constructively and destructively, as implied by Figure 5.
Note that, for zero thickness screens, the predictions by the Vainstein-based model for Im[β 0 ] deviate from those found numerically by about ∼100% for Scale-1  Table II for details.
(M = 0.2), ∼52% for Scale-2 (M = 0.07) and ∼17% for Scale-3 (M = 0.02). Clearly, the agreement is better in the limit of small M (highly-diffractive/overmoded structure, with large Frensel number) where Vainstein's model is expected to increasingly hold. Within this validity limit, one may use the Vainstein-model for a quick first approximation of propagation loss and then add a correction to account for screen thickness effect, as predicted by the presented mode-matching results. For example, at 3 THz, the iris line originally proposed by reference [2] for LCLS (with a = 5.5 cm and b = 30 cm) would give power loss as low as 13.8% per 150 m (to reach the Near Experimental Hall at LCLS) using the Vainstein model and assuming infinitely thin screens (see Figure 3). When the screens are around 2-mm thick, for realistic implementation, we can adjust this value using the mode-matching predictions in Table II as an approximate guide, reducing the final power loss estimate to 9.4%. If mechanically feasible, we can also use a larger iris radius of 10 cm (exploiting the 1/a 3 dependence of loss factor on radius), to achieve an even lower loss of 2.4% per 150 m for infinitely thin screens, which becomes 1.6% for 2-mm thick screens (assuming similar δ effect as in the previous case); see Figure 9.
Using the same approach, Table III summarizes  The influence of screen thickness on field polarization and amplitude profile is shown in Figure (

D. Ohmic losses
Given the small skin depth of good conductors, such as copper or aluminum, at THz frequencies, we can use the well-known perturbation technique to compute the ohmic (conductor) losses at metallic screens, using our knowledge of the tangent magnetic field distribution near the screens when such screens are assumed to be perfectly conductive [14]. Specifically, the attenuation constant due to the conductive loss, α, is calculated by the usual formula where P z = 1 2 R s |H tan | 2 adθdz is the power lost to conductive edges of the screen per unit length of the line, is the total power on the line, and R s = ωµ/(2σ) is the surface resistance for a screen with conductivity σ. Finding the fields through the modematching equations (48)-(59) then substituting in α, we find that the ohmic losses on the iris line are negligible compared to the loss due to diffraction or radiation. For example, the attenuation constant α for a Scale-3 irisline structure with 2-mm thick screens is found to be approximately 2.7 × 10 −7 (1/m) for copper screens and 3.5 × 10 −7 (1/m) for aluminum screens, which are more than 1000-times smaller than the attenuation constant Im[β 0 ] ∼ =0.00041 (1/m) due to diffraction (see Table II). E. Limiting case test (∆ → 0) As a compatibility check, we drove the current modematching model toward the limit of vanishing gap (approaching a smooth line) to compare the predicted eigenvalues with the well-known values for a smooth circular waveguide. As expected, the hybrid mode decoupled into TE and TM modes in this limit, with the values of β 0 now approaching a purely real value (due to suppressed radiation by the closing gaps) and equal to the the analytical values of β TE = k 2 0 − (V mj /a) 2 and β TM = k 2 0 − (V mj /a) 2 , where V mj is the j-th root of the J m function, and V mj is the j-th root of the J m function, [14].
Specifically, using a Scale-1 structure with a gap as small as 2∆ = 10 −5 mm, the mode-matching code gave β TE = 62742.6114531 + i0.0023 and β TM = 62444.425897 + i0.0055, compared to the analytic pre-diction (with at 2∆ = 0) of 62742.611427 + i0 and 62444.425854 + i0, respectively. Notice the good agreement and that the imaginary parts given by the modematching are much smaller (negligible) compared to their values for the open gaps, which had Im[β 0 ] ∼ 26 (see Table II).

F. Transients, mode launching and other considerations
As was mentioned in Section I, the full implementation of this iris-line structure will require addressing further questions, in addition to those discussed above (propagation loss, polarization purity and amplitude profile). These questions pertain to the line dispersion, mode launching and matching from source to line, sensistivity to mechanical tolerances and fabrication misalignments, and the transient regime at the line's entrance (before the mode settles). Such questions are not the subject of this paper and will be examined in a subsequent publication by the authors. Here we briefly remark that mode launching is expected to have a strong influence on the length of the transient at the line's entrance. Indeed, if the incident wave is assumed to be an abstract plane wave, then reaching an equilibrium in the form of a dipole mode on the iris line is expected to be roughly around a distance ∼πa 2 /λ 0 from the entrance (e.g. see [23] and [24]). However, if the field is launched using a specific mode from the source (e.g. a Gaussian profile with linear polarization from a THz linear undulator), then the transient is expected to be much shorter.

V. CONCLUSIONS
Vector field analysis using mode matching and clustered expansions has been presented for an overmoded iris-line structure for THz radiation transportation. The analysis includes the effect of finite screen thickness on the propagation properties of the dominant hybrid mode. As a specific application, the iris line was investigated as a method for THz radiation transport at LCLS, in the frequency range 3-15 THz. It is observed that having finite screen thicknesses, up to one tenth of the period of the iris-line (δ < b/10), seems to reduce the attenuation constant Im[β 0 ] down to approximately −40% for the highly-overmoded THz iris-line structure proposed for THz transport at LCLS. The horizontal polarization purity across the gap seem to be unaffected by the finite screen thickness. The field amplitude profile, which approximately follows a J 0 function profile (close to Gaussian), also seems to be unaffected by the finite screen thickness. The Vainstein-based model predictions, for zero-thickness screens, seems to hold increasingly in agreement with mode-matching numerical predictions for highly-diffractive/overmoded structures (limit of small To nontrivially solve these two equations for the eigenvalue k t (the transverse wavenumber), let us multiply (A8) by E 0 /H 0 and then add/substrate (A9), respectively, to have Let's use perturbation theory in the small paramater M , and expand k t as For the unperturbed case when M → 0, the conditions (A10) and (A11) reduce to From (A13), it is clear that we must either have E 0 = H 0 or E 0 = −H 0 as solutions. The former of these two is the desirable (balanced) hybrid mode that has a peak amplitude profile at the center of the iris and minimum at the walls (akin to the HE11 mode in traditional microwave corrugated waveguides [8,12,13]). The latter of the two, on the other hand, will result in a surface mode that "sticks" to the walls, with a null at the center of the iris, and is undesirably lossy (akin to EH11 mode in traditional microwave corrugated waveguides [8,12,13]). Substituting in (A14) we immediately see how having E 0 = H 0 in (A14) will result in J n−1 (k t a) = 0 with k t0 = V (n−1)j /a, whereas for E 0 = −H 0 we have J n+1 (k t a) = 0 with k t0 = V (n+1)j /a, where V nj denotes the j th root (zero) of the Bessel function J n . The dominant (lowest-order) modes in each of these cases are obtained when n = 1 and give J 0 (2.4r/a) versus J 2 (5.1r/a) profiles, respectively. Note that in the case for n = 0 it can be easily shown that the hybrid mode is degenerated to the TE or TM mode families. Let us now continue our treatment with the desirable balanced mode, with E 0 = H 0 , and for β ≈ k 0 for paraxial high-frequency propagation above cutoff, and using perturbation theory to the first order in M and with k t0 = V (n−1)j /a.
Using the following Bessel identities [21] J n (k t0 a) + J n+1 (k t0 a) = n k t0 a J n (k t0 a), (A15) and exploiting the fact that our boundary condition has already forced J n−1 (k t0 a) = 0 in the unperturbed limit, we can deduce the following relations and perturbative expansions, up to O[M ], Substituting from (A20)-(A22) into the original condition in (A11), to first order in M , and using the shorthand notation J n (k t0 a) ≡ J n and n(n+1) a 2 k 2 t 0 − 1 ≡ ψ, we find c 1 as follows 0 = [J n (k t a) + (1 + i)β 0 aM k t J n (k t a)] + n ak t J n (k t a) where n = 1, 2, 3, · · · .
The lowest-order (dominant) hybrid mode (n = 1, j = 1) will therefore have approximately the desired profile of a J 01 2.4r For small M , this is similar to the usual profile of the J 0 (2.4r/a) function, but without allowing the "skirt" of the J 0 function to go exactly to zero on the radius r = a (as it would have in a smooth pipe). This allows the hybrid mode to exist at the boundary due to coupling (impedance) and gives the desired features of this mode of the iris line. The important result (A23) agrees with that predicted by [2] and is similar to the theory of the HE11 balanced hybrid-mode in traditional microwave corrugated guides (e.g. see [8,12,13]). One should note, however, that the current model using Vainstein's boundary condition is mainly governed by diffraction phenomena, which dominates the propagation loss (A24), rather than by ohmic (conductor) losses.
The propagation loss on the line can now be quantified using the imaginary part of β, as the attenuation constant. Using the fact that β 2 = k 2 0 − k 2 t and the result in (A23), we can now find the power loss (L p ) on the line by simply calculating 1 − e −2Im [β] , which gives, after algebraic manipulation, the following where m ≡ n − 1 = 0, 1, 2, · · · . Equation (A24) indicates lower losses for higher frequencies (as a function of ω −1/2 ) and for larger radii (as a function of a −3 ). Equation (2) is a special case of (A24), taken at n = 1, j = 1.
We can now proceed to examine the polarization of this desirable mode. For this balanced mode with n = 0, we find that (A2) can be reduced to To examine the E field polarization, let us now convert equations (A25) and (A26) to the Cartesian coordinates using the transformationr =x cos θ +ŷ sin θ and θ = −x sin θ +ŷ cos θ. After simplifying, this gives the following expression for the transverse field This result highlights the character of the polarization for the iris-line balanced modes. For the dominant mode of interest (n = 1, j = 1), (A29) and (A23) give which means that the field is only polarized horizontally and it remains so across the entire aperture of the iris (fixed dipole polarization) while its amplitude profile is approximately ≈ J 0 (2.4r/a). This gives an intensity profile that is similar to a Gaussian intensity (note: the profile of J 2 0 has about 97% overlap with a perfect Gaussian intensity, when both are normalized to the same maximum and have same 1/e 2 width below the maximum; see Figure 11b). Therefore, the dominant balanced mode of the iris-line is suitable for direct coupling with undulator radiation (see Figure 11a). As noted earlier, the complex argument of the Bessel function, with small complex perturbation proportional to M will cause the edges of the function skirts to be slightly lifted off zero at r = a, allowing for the hybrid nature of the mode to be established. Indeed, this can be seen explicitly by writing out the real part of the final field expression. Since the Bessel's function J n is an entire function (i.e. analytic everywhere in the complex plane), the profile of E x11 in (A31) can be explicitly reduced further (with M > 0) to its real and imaginary parts. For any fix r value, one can Taylor expand in the small parameter M to find that the field profile is proportional to Thus, the real part of field profile approximately given as which is shown in Figure 12. This shows the polarization and amplitude profile features of the dipole mode on the iris line. As one may expect, this mode is akin to dipole modes traditionally seen in similar periodic structures when driven by ultrarelativistic electron bunchs; e.g. see [15,25].
Appendix B: Formal justification for the use of clustering to approximate field expansions in an overmoded paraxial line In this appendix we use perturbation theory to derive a justification to the practice of clustering in field expansions, as used in Section III for a highly-overmoded iris line with paraxial incidence. Specifically, we wish to show that one can assume that the following two summations are interchangeable where the first sum on the RHS is clustered around 0 (dominant harmonic, at n=0), while the second sum is clustered around its image (at n=−n 0 =−2N 0 ), where N 0 is defined as the nearest integer to the number of wavelengths per structure period, b/λ 0 ; see Figure 13. Here nSteps is chosen by us to provide a local (trunacted) span around the dominant harmonic and its image, such that the overall summation is shorter than the regular summation that would run continuously from −N max to +N max , with the arbitrary index N max chosen high enough as to encompass all harmonics between and around the dominant harmonic and its image. The idea behind such a practice stems from an intuitive assumption that the dominant mode corresponds to an almost-TEM wave that is paraxially incident onto the overmoded iris line and slightly perturbed by the presence of the large irises, as discussed in Section III. Let's show that this assumption is formally justifiable. Without loss of generality, we specialize the discussion to the dipole mode (m = 1) of the iris line. We know that the propagation constant for the n th harmonic is β n = β 0 + 2πn/b, where β 0 approaches the wavenumber k 0 for paraxial propagation high above cutoff. In the limit of an ideal line (lossless with TEM wave), we can put β n = k 0 + 2πn/b, which would lead to a transverse wavenumber k tn = k 2 0 − β 2 n = i (2πn/n) 2 + 4πnk 0 /b. The dipole mode field equations at a given point along the line, say z = 0, are given [from It is clear from k tn = i (2πn/n) 2 + 4πnk 0 /b that certain n values will lead to real k tn values, while other n values will lead to purely imaginary k tn = iγ n , where γ n is real. For the former, the Bessel functions in (B2)-(B7) are Bessel functions of the first type [e.g. J 1 (k tn r)], whereas for the latter case, they become, using the relation J m (iγ n r) = i m I m (γ n r), modified Bessel functions of the first type [e.g. I 1 (γ n r)]. To work in terms of real k tn we will need negative n values that satisfy 4π|n|k 0 /b ≥ 4|n| 2 π 2 /b 2 , which leads to the requirement that −n 0 ≤ n ≤ 0, where n 0 = 2b/λ 0 = 2N 0 . For the remaining values of n, γ n is real. See Figure 13. Given the symmetry between the two types of formulations, we can work in either set of wavenumbers and Bessel functions (k tn and J m , or γ n and I m ) and expect the other to have similar features. We can therefore confine the remainder of this discussion to the formulation using k tn and J m functions over the middle band shown in Figure 13. For the same iris-line dimensions used in Section IV, we have a = 0.05 m, b = 0.33 m, λ 0 = 10 −4 m, which give −6666 = −n 0 ≤ n ≤ 0. We recognize the edges of this range as the terms we initially assumed to be the dominant harmonic and its image. Let's next show that the strength of all the harmonics between these two terms (n = 0, −n 0 ) get smaller and smaller as we move away from edges.
Consider the total complex power associated with each harmonic travelling down the iris line, by calculating the integral (call it I) of the Poynting vector over the line's cross section, where the asterisk denotes complex conjugation. Before substituting the modal field expressions (B2)-(B7) into (B8) to find the respective powers carried by the harmonics, we should lift the idealistic assumption of lossless line, since the line will exhibit some loss. We introduce two perturbation parameters and write the complex propagation constant as where ε 1 1 and ε 2 1 represent, respectively, the small shift in real part of β 0 relative to k 0 and the small attenuation experienced by the waves on the line. Working to first order in the small parameters ε 1 , ε 2 , this leads, after algebraic manipulation, to Within the present range of n, we define the convenient positive-integer measure d as d = n 0 − |n|, using it as our variable, and we substitute n → −|n| in (B10)-(B13), to yield, after reduction to first order in ε 1 , ε 2 and renaming φ n as φ d ,