The least square progressive iterative approximation (LSPIA) is an efficient method for large-scale data fitting. To accelerate the convergence rate of LSPIA, we propose two kinds of LSPIA for curves and surfaces based on Newton iterative method. Firstly, we construct a multivariate function with control points as variables, whose Hessian matrix is a positive definite matrix. The control points corresponding to the minimal value of the multivariate function are consistent with the convergence results of LSPIA. Meanwhile, the Newton iterative method is applied to solve the minimal value of the multivariate function. Then we obtain two LSPIA iterative schemes by using singular value decomposition and the conjugate gradient method respectively to deal with the Hessian matrix and adjustment vectors, which are named NLSPIA and INLSPIA respectively. Finally, we also give the theoretical proof of the convergence of the proposed schemes. Numerical examples verify the effectiveness and feasibility of the proposed methods, and also show that the convergence rate and the calculation time of NLSPIA and INLSPIA are both better than those of the classical LSPIA under the same conditions.