  We present an algorithm for calculating the metric perturbations and
gravitational self-force for extreme-mass-ratio inspirals (EMRIs) with
eccentric orbits. The massive black hole is taken to be Schwarzschild and
metric perturbations are computed in Lorenz gauge. The perturbation equations
are solved as coupled systems of ordinary differential equations in the
frequency domain. Accurate local behavior of the metric is attained through use
of the method of extended homogeneous solutions and mode-sum regularization is
used to find the self-force. We focus on calculating the self-force with
sufficient accuracy to ensure its error contributions to the phase in a long
term orbital evolution will be $\delta\Phi \lesssim 10^{-2}$ radians. This
requires the orbit-averaged force to have fractional errors $\lesssim 10^{-8}$
and the oscillatory part of the self-force to have errors $\lesssim 10^{-3}$ (a
level frequently easily exceeded). Our code meets this error requirement in the
oscillatory part, extending the reach to EMRIs with eccentricities of $e
\lesssim 0.8$, if augmented by use of fluxes for the orbit-averaged force, or
to eccentricities of $e \lesssim 0.5$ when used as a stand-alone code. Further,
we demonstrate accurate calculations up to orbital separations of $a \simeq 100
M$, beyond that required for EMRI models and useful for comparison with
post-Newtonian theory. Our principal developments include (1) use of fully
constrained field equations, (2) discovery of analytic solutions for
even-parity static modes, (3) finding a pre-conditioning technique for outer
homogeneous solutions, (4) adaptive use of quad-precision and (5) jump
conditions to handle near-static modes, and (6) a hybrid scheme for high
eccentricities.
