The $f$-divergence is a fundamental notion that measures the difference between two distributions. In this paper, we study the problem of approximating the $f$-divergence between two Ising models, which is a generalization of recent work on approximating the TV-distance. Given two Ising models $\nu$ and $\mu$, which are specified by their interaction matrices and external fields, the problem is to approximate the $f$-divergence $D_f(\nu\,\|\,\mu)$ within an arbitrary relative error $\mathrm{e}^{\pm \varepsilon}$. For $\chi^\alpha$-divergence with a constant integer $\alpha$, we establish both algorithmic and hardness results. The algorithm works in a parameter regime that matches the hardness result. Our algorithm can be extended to other $f$-divergences such as $\alpha$-divergence, Kullback-Leibler divergence, R\'enyi divergence, Jensen-Shannon divergence, and squared Hellinger distance.