Ising model density of states 1D Ising model: exact result
In the 1D Ising model with fixed $J_{ij} = J$, without magnetic field, the density of states (dos) can be calculated exactly. There is a caveat in the case of periodic boundaries (which I don't expect to change the result considerably), but let us consider a linear spin chain. 
The Hamiltonian (taking $J=1$)
$$
H = - \sum\limits_{\langle i,j \rangle} \sigma_i \sigma_j
$$
takes values $-N+1, -N+3, \ldots, N-3, N-1$ for a system with $N$ spins (and periodic boundaries). For $z=2$ this has the special property that the number of terms in the Hamiltonian is $\frac{(N-1)z}{2} = N-1$, one less than the number of spins because the two ends have no interaction. Note that the terms $H_{i,i+1} \equiv -\sigma_i \sigma_{i+1}$ are linearly independent. In fact,  there is a 2-to-1 correspondence between the set of spins $\{\sigma_i\}$ and the set of interaction terms $\{H_{ij}\}$. Each set of interactions is invariant under $\sigma_i \rightarrow -\sigma_i$ (for all spins). Therefore, the number of states with given energy $E(k) = -N+1+2k$ with $k \in \{0, 1, \ldots, N-1\}$ is given by
$$
g(E(k)) \equiv g(k) = 2 \binom{N-1}{k}
$$
Indeed, summing over all $k$ gives
$$
\text{Total number of states} = \sum\limits_{k=0}^{N-1} 2 \binom{N-1}{k} = 2^N.
$$
Higher dimensions: are there analytical approximations?
Now I've been looking into ways to extend such an argument to higher dimensions, but so far found nothing. The problem is that for $d>1$ there is no neat correspondence between the spins and the interaction terms. For instance, given two spins on a 2D lattice, there are multiple paths between the two, so the interaction terms are limited by the constraint that the spin values are consistent over each of those paths. If one path predicts spin 1 should have an opposite sign from spin 2 and the other that the signs should be equal, then such a set of interactions is not permitted. 
I know that there is a neat algorithm called the Wang-Landau algorithm for calculating the density of states. In fact, the results for the 2D Ising model look similar to the binomial distribution, though appearing somewhat steeper (see the original article). However, I was not able to find any analytical work on calculating the dos. Are there any methods that can be used to obtain an (approximate) analytical result?
 A: If you look for an exact density of state (DOS) for a 2D Ising model I recommend the Paul Beale's article PhysRevLett.76.78.
In the article the author shows how to correctly obtain the density of states starting from the definition of the partition function
$Z_{n,m}(K)=e^{2nmK}\sum_{k=0}^{nm}g_kx^{2k}$,
where n is the number of columns, m is the number of
rows, $K=J/k_BT$ is the coupling, $k_B$ is Boltzmann’s
constant, $x=e^{-2K}$ is the low-temperature expansion
variable, and the coefficient $g_k$ is the number of config-
urations with energy $4kJ$ above the two ground states (all
spins +1 or all spins -1).
The result of the partition function has a binomial term that represents the density of state. The solution has many other terms and your calculation is very hard work. So the author give a link for download a notebook from Wolfram Matematica with the implementation of the calculations of the density of states and other thermodynamics quantities.
In the article the link is broken, but I found another in the author’s personal website. I hope to have help with this question.
Link for download the notebook
