  The past decade has seen significant progress in understanding galaxy
formation and evolution using large-scale cosmological simulations. While these
simulations produce galaxies in overall good agreement with observations, they
employ different sub-grid models for galaxies and supermassive black holes
(BHs). We investigate the impact of the sub-grid models on the BH mass
properties of the Illustris, TNG100, TNG300, Horizon-AGN, EAGLE, and SIMBA
simulations, focusing on the M_BH-M_star relation and the BH mass function. All
simulations predict tight M_BH-M_star relations, and struggle to produce the
lowest (M_BH< 10^7.5 Msun) in galaxies of M_star~10^10.5-10^11.5 Msun. While
the time evolution of the mean M_BH-M_star relation is mild (<1 dex in BH mass
for 0<z<5) for all the simulations, its linearity (shape) and normalization
varies from simulation to simulation. The strength of SN feedback has a large
impact on the linearity and time evolution for M_star<10^10.5 Msun. We find
that the low-mass end is a good discriminant of the simulation models, and
highlights the need for new observational constraints. At the high-mass end,
strong AGN feedback can suppress the time evolution of the relation
normalization. Compared with the observations of the local universe, we find an
excess of BHs with M_BH>10^9 Msun in most of the simulations. The BH mass
function is dominated by efficiently accreting BHs (log10 f_Edd >-2$) at high
redshifts, and transitions progressively from the high-mass to the low-mass end
to be governed by inactive BHs. The transition time and the contribution of
active BHs are different among the simulations, and can be used to evaluate
models against observations.
