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 less than or similar to 10(-2) rad. This requires the orbit-averaged force to have fractional errors less than or similar to 10(-8) and the oscillatory part of the self-force to have errors less than or similar to 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 less than or similar to 0.8, if augmented by use of fluxes for the orbit-averaged force, or to eccentricities of e less than or similar to 0.5 when used as a stand-alone code. Further, we demonstrate accurate calculations up to orbital separations of a similar or equal to 100M, 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 preconditioning technique for outer homogeneous solutions, (4) adaptive use of quad precision, (5) jump conditions to handle near-static modes, and (6) a hybrid scheme for high eccentricities.
