Multi-fidelity surrogate model-based engineering optimization has received much attention because it alleviates the computational burdens of expensive simulations or experiments. However, due to the nonlinearity of practical engineering problems, the initial sample set selected to produce the first set of data will almost inevitably miss certain features of the landscape, and thus, the construction of a useful surrogate often requires further, judicious infilling of some new samples. Sequential sampling strategies used to select new infilling samples during each iteration can gradually extend the data set and improve the accuracy of the initial model with an acceptable cost. In this paper, a sequential sampling generation method based on the Voronoi region and the sample density (SSGM-VRDS) is proposed. First, with a Monte Carlo-based approximation of a Voronoi tessellation for region division, Pearson correlation coefficients and cross-validation (CV) are employed to determine the candidate Voronoi region for infilling a new sample. Then, a relative sample density is defined to identify the position of the new infilling point at which the sample is the sparsest within the selected Voronoi region. A correction of this density is carried out concurrently through an expansion coefficient. The proposed method is applied to three numerical functions and a lightweight design problem via finite element analysis (FEA). Results suggest that the SSGM-VRDS strategy has outstanding effectiveness and efficiency in selecting a new sample for improving the accuracy of a surrogate model, as well as practicality for solving practical optimization problems.