We investigate the abundance of galactic molecular hydrogen (H$_2$) in the Evolution and Assembly of GaLaxies and their Environments (EAGLE) cosmological hydrodynamic simulations. We assign H$_2$ masses to gas particles in the simulations in post-processing using two different prescriptions that depend on the local dust-to-gas ratio and the interstellar radiation field. Both result in H$_2$ galaxy mass functions that agree well with observations in the local and high-redshift Universe. The simulations reproduce the observed scaling relations between the mass of H$_2$ and the stellar mass, star formation rate and stellar surface density. Towards high edshifts, galaxies in the simulations display larger H$_2$ mass fractions, and correspondingly lower H$_2$ depletion timescales, also in good agreement with observations. The comoving mass density of H$_2$ in units of the critical density, $Omega_{rm H_2}$, peaks at $zapprox 1.2-1.5$, later than the predicted peak of the cosmic star formation rate activity, at $zapprox 2$. This difference stems from the decrease in gas metallicity and increase in interstellar radiation field with redshift, both of which hamper H$_2$ formation. We find that the cosmic H$_2$ budget is dominated by galaxies with $M_{rm H_2}>10^9,rm M_{odot}$, star formation rates $>10,rm M_{odot},rm yr^{-1}$ and stellar masses $M_{rm stellar}>10^{10},rm M_{odot}$, which are readily observable in the optical and near-IR. The match between the H$_2$ properties of galaxies that emerge in the simulations and observations is remarkable, particularly since H$_2$ observations were not used to adjust parameters in EAGLE.
Download