Efficient and accurate numerical simulation of seismic wave propagation is important in various Geophysical applications such as seismic full waveform inversion (FWI) problem. However, due to the large size of the physical domain and requirement on low numerical dispersion, many existing numerical methods are inefficient for numerical modelling of seismic wave propagation in an heterogenous media. Despite the great efforts that have been devoted during the past decades, it still remains a challenging task in the development of efficient and accurate finite difference method for multi-dimensional acoustic wave equation with variable velocity. In this paper we proposed a Padé approximation based finite difference scheme for solving the acoustic wave equation in threedimensional heterogeneous media. The new method is obtained by combining the Padé approximation and a novel algebraic manipulation. The efficiency of the new algorithm is further improved through the Alternative Directional Implicit (ADI) method. The stability of the new algorithm has been theoretically proved by energy method. The new method is conditionally stable with a better Courant -Friedrichs -Lewy condition (CFL) condition, which has been verified numerically. Extensive numerical examples have been solved, which demonstrated that the new method is accurate, efficient and stable.