A non-vanishing neutrino mass leads to a few percent suppression of matter perturbations, with respect to a massless neutrino cosmology, on all scales below their free-streaming scale. The effect allows to place constraints of the sum of neutrino masses from measurements of the galaxy power spectrum but it also constitutes a theoretical uncertainty in dark energy studies. For both these reasons its accurate description is particularly relevant. I will review how well numerical simulations and analytical methods can account for neutrino masses and argue that the interplay between the two approaches is particularly significant and fruitful for this specific problem.