A perfectly electrically conducting (PEC) empty waveguide, of constant cross-section and axis aligned with the z-axis, supports time-harmonic propagating modes of the form u(x,y)e −1ωt–iγz , where the function u satisfies the two-dimensional Helmholtz equation (Δ + λ)u = 0, k denoting wavenumber and λ=k 2 –γ 2 . Propagating modes occur at values of λ where non-trivial solutions of the Helmholtz equation occur; the cutoff wavenumber for each mode corresponds to setting γ to zero. If axially aligned PEC structures are inserted in the empty waveguide, the field distribution of each propagating mode is distorted and the corresponding propagation constant λ i is perturbed to a value λ i + △λ i . Some estimates of the perturbation △λ i obtained for particular structures are given in [1]. Accurate calculation of the perturbation to values of the propagation constants presents significant numerical challenges, particularly if the insert has sharp edges or comers, or is touching the inner boundary of the waveguide (or nearly so). Equally, accurate calculation of the modal field distribution is required in order to track its evolution from the empty waveguide mode as the size of the insert (or inserts) increases. In this context, any satisfactory method must provide cutoff wave numbers with any pre-specified accuracy and must be capable of resolving distinct modes in the spectrum even when their cutoff frequencies lie extremely close.