We solve a singe-particle Dirac equation with Woods-Saxon potentials using an iterative method in the coordinate space representation. By maximizing the expectation value of the inverse of the Dirac Hamiltonian, this method avoids the variational collapse, in which an iterative solution dives into the Dirac sea. We demonstrate that this method works efficiently, reproducing the exact solutions of the Dirac equation.