Gaussian Mixture Models (GMMs) are important tools for modeling complex data in machine learning tasks and computer vision applications. However, computing f-divergences between GMMs remains challenging due to the absence of a closed-form expression, which led to expensive numerical approximations that limit practical applications. In this paper, we give an efficient f-divergence approximation through the embedding of GMMs into the symmetric positive definite (SPD) matrices. Our main result is that for any compact set of non-degenerate GMM parameters, the f-divergence between two GMMs and the computationally efficient f-divergence between their corresponding centered multivariate normal distributions in the SPD space are uniformly equivalent. Our approach preserves the geometric structure of GMMs while enabling closed-form computation. As an instance the proposed framework is applied on the UIUC texture recognition datasets.