  We present direct logarithmically optimal in theory and fast in practice
algorithms to implement the tensor product high order finite element method on
multi-dimensional rectangular parallelepipeds for solving PDEs of the Poisson
kind. They are based on the well-known Fourier approaches. The key new points
are the fast direct and inverse FFT-based algorithms for expansion in
eigenvectors of the 1D eigenvalue problems for the high order FEM. The
algorithms can further be used for numerous applications, in particular, to
implement the tensor product high order finite element methods for various
time-dependent PDEs. Results of numerical experiments in 2D and 3D cases are
presented.
