The design of parallel algorithms and architectures for solving linear systems using two-step division-free Gaussian elimination method is considered. The two-step method circumvents the ordinary single-step division-free method by its greater numerical stability. In spite of the rather complicated computations needed at each iteration of the two-step method, we develop first an innovative regular iterative algorithm, then a two-dimensional array processor by deriving a localized dependency graph of the algorithm and adopting a systematic approach to investigate the set of all admissible solutions and obtain the optimal architecture under linear scheduling. The optimal array processor improves the previous systolic designs based on the widely used Gaussian elimination in term of numerical stability and the time-space complexity for VLSI implementation because of the absence of division operations.