Molecular hydrogen ($ H_2 $) plays a crucial role in the formation and evolution of galaxies, serving as the primary fuel reservoir for star formation. In a metal-enriched Universe, $ H_2 $ forms mostly through catalysis on interstellar dust grain surfaces. However, due to the complexities of modelling this process, star formation in cosmological simulations often relies on empirical or theoretical frameworks that have only been validated in the local Universe to estimate the abundance of $ H_2 The goal of this work is to model the connection between the processes of star, dust, and $ H_2 $ formation in our cosmological simulations. Building upon our recent integration of a dust evolution model into the star formation and feedback model MUPPI, we included the formation of molecular hydrogen on the surfaces of dust grains. We also accounted for the destruction of molecules and their shielding from harmful radiation. The model reproduces, reasonably well, the main statistical properties of the observed galaxy population for the stellar, dust, and components. The evolution of the molecular hydrogen cosmic density ($ H2 $) in our simulated boxes peaks around redshift $z=1.5$, consistent with observations. Following its peak, $ H2 $ decreases by a factor of two towards $z=0$, which is a milder evolution than observed. Similarly, the evolution of the molecular hydrogen mass function since $z=2$ displays a gentler evolution when compared to observations. Our model recovers satisfactorily the integrated molecular Kennicut-Schmidt (mKS) law between the surface star formation rate ($ SFR $) and surface density ($ H2 $) at $z=0$. This relationship is already evident at $z=2$, albeit with a higher normalization. We find hints of a broken power law with a steeper slope at higher $ H2 $. We also study the $ H_2 $-to-dust mass ratio in galaxies as a function of their gas metallicity and stellar mass, observing a decreasing trend with respect to both quantities. The $ H_2 $-to-dust mass fraction for the global population of galaxies is higher at higher redshift. The analysis of the atomic-to-molecular transition on a particle-by-particle basis suggests that gas metallicity cannot reliably substitute the dust-to-gas ratio in models attempting to simulate dust-promoted $ H_2