Abstract The standard deviations of the observation function can be calculated by the law of covariance propagation. However, the approximate function method based on Taylor series expansion requires complex derivative operations when solving the standard deviation of nonlinear function. The Monte Carlo method can avoid the derivative operation, but it is not objective in the selection of the number of simulations. Moreover, the Monte Carlo method cannot directly control the simulation results. To overcome these disadvantages, we introduce the Stein two-stage method into the covariance propagation theory of nonlinear functions. Combined with the Monte Carlo method, we design the Stein Monte Carlo (SMC) algorithm flow of nonlinear function covariance propagation. We apply the SMC method in the two-dimensional polynomial function and covariance propagation of GNSS baseline vector. Results verify the effectiveness of the SMC method. This method provides a new idea for the covariance propagation of nonlinear models.