A high-precision three-dimensional forward simulation method based on optimized finite-difference coefficients is proposed to address the numerical dispersion problem in three-dimensional finite-difference forward simulation. Firstly, the spatio-temporal dispersion relationship is derived from the three-dimensional acoustic wave equation, and an infinite norm optimization problem is constructed using dispersion errors. Then, a new objective function is obtained by utilizing the symmetry of the error function. Finally, the optimization finite-difference coefficients which will vary with velocity is obtained by solving the objective function through Remez iteration method. Dispersion analysis shows that the optimized finite-difference coefficients obtained by the new method have the characteristic of equal ripple error, which means the error in the low wave number region changes regularly within the error limit. Compared with the finite-difference coefficients obtained in the spatial domain, the finite-difference coefficients obtained by the new method can effectively adapt to changes in the CFL (Courant Friedrichs Lewy) number, and can suppress time dispersion under large CFL number. Compared with the finite-difference coefficients obtained by the Taylor expansion method in the spatial domain, the finite-difference coefficients obtained by the new method have a wider frequency band and can reduce the spatial dispersion caused by the errors in high wave number region. The high-speed homogeneous model example shows that the time dispersion of the new method is invisible and the spatial dispersion is slight in high-speed formations, proving that the new method can suppress both time dispersion and spatial dispersion and has a high accuracy in simulating wave fields. The salt model shows that the new method can significantly suppress numerical dispersion caused by both low-speed layers and high-speed geological bodies in strata with rapidly changing velocity, significantly improve the resolution of seismic events in seismic records, and verify that the new method is suitable for complex models with drastic velocity changes.