  In this work we present RESCU, a powerful MATLAB-based Kohn-Sham density
functional theory (KS-DFT) solver. We demonstrate that RESCU can compute the
electronic structure properties of systems comprising many thousands of atoms
using modest computer resources, e.g. 16 to 256 cores. Its computational
efficiency is achieved from exploiting four routes. First, we use numerical
atomic orbital (NAO) techniques to efficiently generate a good quality initial
subspace which is crucially required by Chebyshev filtering methods. Second, we
exploit the fact that only a subspace spanning the occupied Kohn-Sham states is
required, and solving accurately the KS equation using eigensolvers can
generally be avoided. Third, by judiciously analyzing and optimizing various
parts of the procedure in RESCU, we delay the $O(N^3)$ scaling to large $N$,
and our tests show that RESCU scales consistently as $O(N^{2.3})$ from a few
hundred atoms to more than 5,000 atoms when using a real space grid
discretization. The scaling is better or comparable in a NAO basis up to the
14,000 atoms level. Fourth, we exploit various numerical algorithms and, in
particular, we introduce a partial Rayleigh-Ritz algorithm to achieve
efficiency gains for systems comprising more than 10,000 electrons. We
demonstrate the power of RESCU in solving KS-DFT problems using many examples
running on 16, 64 and/or 256 cores: a 5,832 Si atoms supercell; a 8,788 Al
atoms supercell; a 5,324 Cu atoms supercell and a small DNA molecule submerged
in 1,713 water molecules for a total 5,399 atoms. The KS-DFT is entirely
converged in a few hours in all cases. Our results suggest that the RESCU
method has reached a milestone of solving thousands of atoms by KS-DFT on a
modest computer cluster.
