The phonon density of states tells us how many phonon modes there are at every frequency. This can be calculated by choosing a uniform grid of $\vec{k}$ states in the first Brillouin zone and calculating the frequencies of the phonon modes for each $\vec{k}$ state. A histogram is then constructed out of these frequencies to plot the density of states.

For a bcc lattice, the Brillouin zone is a rhombic dodecahedron. The rhombic dodecahedron can be divided into 48 identical polyhedra that have corners at $\Gamma$, $N$, $H$, and $P$.

$\frac{C_2}{C_1} = $  

$D(\omega)$
	$\large \sqrt{\frac{M}{C_1}}\omega$

There is still some numerical noise visible in the density of states because the $\vec{k}$ state grid is not fine enough. To get a smooth density of states, it would be necessary to calculate more points than is possible in the time it takes to load a web page. Notice that the number of numerical errors is listed. The code that calculates the eigenvalues of the matrix sometimes fails. This is common for simple routines that calculate eigenvalues. As long as there are not too many errors compared to the number of $\vec{k}$ states calculated, it should not have a large influence on the result.

Matlab code for performing this calculation more precisely can be found here.