The linear error propagation and its system matrix of the spatial discretization of matrix SPH, which was derived through perturbation method, were used to analyze the stability of SPH. Based on the linear form, the sufficient condition, which was ever obtained by Swegle, was also derived. Omitting the effect from the continuity equation and analyzing the characteristic equation of the system matrix, two matrixes, which having equivalent eigenvalues to those of the system matrix, were obtained and could respectively represent the existences of the tensile instability and the high frequency instability (HFI) in SPH. Based on the matrix that represents the existence of tensile instability, the tensile instability was found to be sensitive to the phase differences among the errors on particles. Based on other matrix that represents the existence of HFI of SPH, the origin of HFI, i.e. the rank deficiency in stiff matrix, was further found to be correlated with the asymmetry of the stiff matrix induced by the heterogeneous particle volumes. At last, another smoothing length refreshing model was proposed and was tested in two numerical cases.
