  We present a quasi-linearly scaling, first order polynomial finite element
method for the solution of the magnetostatic open boundary problem by splitting
the magnetic scalar potential. The potential is determined by solving a
Dirichlet problem and evaluation of the single layer potential by a fast
approximation technique based on Fourier approximation of the kernel function.
The latter approximation leads to a generalization of the well-known
convolution theorem used in finite difference methods. We address it by a
non-uniform FFT approach. Overall, our method scales O(M + N + N log N) for N
nodes and M surface triangles. We confirm our approach by several numerical
tests.
