Abstract2‐D air‐gap magnetic field distribution is the essential prerequisite of electrical machine analysis. Because of the natural periodicity of rotary machines, Fourier analysis is a suitable choice for air‐gap field prediction. However, due to the end‐effect, periodicity is not present in linear machines and the Fourier series is not appropriate. In engineering mathematics, a rigorous method of solving partial differential equations is the separation of the variables method (SVM) which is based on the hypothesis that the field solution is in the form of the product of two functions of orthogonal directions; for example, in an axisymmetric structure, a longitudinal harmonic function (LHF) and a radial harmonic function (RHF). A particular case of these functions is trigonometric functions, which result in the Fourier series. SVM is imposed in a different manner, where RHF is approximated by the Bessel–Fourier series and consequently LHF has a (piecewise) exponential behaviour. This choice of basis functions not only serves the purpose of modelling the end‐effect but also removes the Gibbs phenomenon, that is, it precisely models discontinuities of flux density at the surface of PMs. A numerical case‐study of a slotless double‐sided linear tubular surface‐PM machine shows that the piecewise exponential approximation can attain 1% error while utilising only three harmonics, whereas the trigonometric approximation, due to the Gibbs phenomenon, slowly arrives at 20% error with 100 harmonics for approximating flux density. Results of both methods were validated using the finite element method (FEM). Additionally, it was discovered that the computational complexity of the model is independent of varying position; therefore, the analytical method could yield the back electromotive force and electromagnetic thrust force for a number of 320 positions at merely 0.2 s while the FEM software required 5 min.