when considering the computation of multipole expansions in particle accelerator magnets, often one has to consider boundary conditions that are accurately represented by a combination of a small number of Fourier modes. Therefore, if a Fourier spectral approach is employed in the azimuthal direction, the resulting computational cost can be significantly reduced.

Concerning specifically particle accelerator simulation, it has been recently shown in [@abele:2019] that the accuracy of particle tracking methods employed to design new accelerator configurations is essentially limited by the accuracy of the field reconstruction. In particular, numerical approximations of electric and magnetic fields that do not respect the inherent mathematical structure of the Maxwell equations may disrupt the effective accuracy of the simulation. This motivates the search for consistent numerical methods that provide high order accuracy in the spatial approximation.

In this work, we propose a combination of a spectral Fourier approximation in the azimuthal direction with an approach in the radial and axial directions, extending the strategy presented in [@oh:2015]. The combination of the Fourier basis and is in the spirit of the recent developments in [@perotto:2017] in computational fluid dynamics. has become popular in applied electromagnetism in the last decade [@bontinck:2017], [@buffa:2010]. One of its 