Abstract: Considering that the conventional geometric iteration method(GIM) has only 1-order convergence, this paper proposes an energy function with second-order smoothness that characterizes the difference between the up-to-date curve and the given data points. In numerical implementation, we initialize the spline curve according to the initial control points and the associated basis functions, then compute the difference function as well as the gradients with regard to the moveable control points, and finally use the L-BFGS technique to find the optimal interpolation/approximation curve. Experimental results show that our accelerated algorithm has a super-linear rate of convergence. With the same accuracy requirement, the improved GIM outperforms the original version by tens to hundreds of times in terms of efficiency. It can be used to both the interpolation problem and the approximation problem, and even to the case where the input data points have variable parameters.