Abstract: As for the ill-posed coefficient matrix in the process of co-seismic slip distribution inversion, a two-step solution method is proposed, and the Laplace second-order smoothing matrix is used for smoothing constraint in this paper. This method not only solves the ill-posed problem of the coefficient matrix, but also suppresses the large gradient variation between adjacent fault patches. During the inversion process of two-step solution method, the L curve method is used to determine the regularization parameters. In this paper, systematic simulation experiments are carried out and the Lushan earthquake is used as an actual earthquake case. The maximum slip inversion results of simulation experiments show that the inversion result of two-step solution method has been improved between 3.34% and 19% than that of one-step solution method in accuracy. The root mean square error inversion result of two-step method is less between 3.3% and 13.3% than that of one-step solution method. The inversion results of the Lushan earthquake show that it is feasible to use the two-step solution for slip distribution inversion.