The two-stage multi-splitting iterative method, established in and developed since 1990, is a kind of new and efficient parallel algorithm for solving the large sparse linear system of equations. As an illustrative example, the new algorithm is applied to the load flow solution of large-scale distribution systems. Through the inner/outer splitting of the Jacobian matrix, the algorithm derived is used for load flow solution of both radial distribution systems and distribution networks with loops. Tests on a 32-node system and a practical distribution system show that the method proposed is robust and efficient for load flow solution of distribution systems. In the stage of the single-splitting of Jacobian matrix now, it has broad prospects for applying the two-stage multi-splitting iterative method in large-scale power system analysis.