Frequency-domain algorithm for the Lorenz-gauge gravitational self-force

State-of-the-art computations of the gravitational self-force (GSF) on massive particles in black hole spacetimes involve numerical evolution of the metric perturbation equations in the time-domain, which is computationally very costly. We present here a new strategy, based on a frequency-domain treatment of the perturbation equations, which offers considerable computational saving. The essential ingredients of our method are (i) a Fourier-harmonic decomposition of the Lorenz-gauge metric perturbation equations and a numerical solution of the resulting coupled set of ordinary equations with suitable boundary conditions; (ii) a generalized version of the method of extended homogeneous solutions [Phys. Rev. D {\bf 78}, 084021 (2008)] used to circumvent the Gibbs phenomenon that would otherwise hamper the convergence of the Fourier mode-sum at the particle's location; and (iii) standard mode-sum regularization, which finally yields the physical GSF as a sum over regularized modal contributions. We present a working code that implements this strategy to calculate the Lorenz-gauge GSF along eccentric geodesic orbits around a Schwarzschild black hole. The code is far more efficient than existing time-domain methods; the gain in computation speed (at a given precision) is about an order of magnitude at an eccentricity of 0.2, and up to three orders of magnitude for circular or nearly circular orbits. This increased efficiency was crucial in enabling the recently reported calculation of the long-term orbital evolution of an extreme mass ratio inspiral [Phys. Rev. D {\bf 85}, 061501(R) (2012)]. Here we provide full technical details of our method to complement the above report.