We transform the seismic wave equations in 2D inhomogeneous anisotropic media into a system of first-order partial differential equations with respect to time t. Based on the transformed equations, a new nearly analytic discrete method (NADM) is developed in this article. Our method enables wave propagation to be simulated in two dimensions through generally anisotropic and heterogeneous models. The space derivatives are calculated by using an interpolation approximation, while the time derivatives are replaced by a truncated Taylor expansion. Our analyses show that the error of the NADM is less than that of the conventional finite-difference method (FDM) and is about 1/60 to 1/100 of that of the FDM. We also demonstrate numerically that the stability of the NADM is higher than that of the FDM. The three-component seismic wave fields in two layered isotropic and transversely isotropic media (TIM) are simulated and compared with the conventional FDM. Again, we show from the three-component seismic wave fields that the NADM has higher accuracy, stronger stability, and less numerical dispersion, effectively suppressing the source noises as compared with the FDM.