We present a simulation of the cosmic evolution of the atomic and molecular phases of the cold hydrogen gas in about 3 X 107 galaxies, obtained by postprocessing the virtual galaxy catalog produced by De Lucia & Blaizot on the Millennium Simulation of cosmic structure. Our method uses a set of physical prescriptions to assign neutral atomic hydrogen (H I) and molecular hydrogen (H2) to galaxies, based on their total cold gas masses and a few additional galaxy properties. These prescriptions are specially designed for large cosmological simulations, where, given current computational limitations, individual galaxies can only be represented by simplistic model objects with a few global properties. Our recipes allow us to (1) split total cold gas masses between H I, H2, and helium, (2) assign realistic sizes to both the H I and H2 disks, and (3) evaluate the corresponding velocity profiles and shapes of the characteristic radio emission lines. The results presented in this paper include the local H I and H2 mass functions, the CO luminosity function, the cold gas mass-diameter relation, and the Tully-Fisher relation (TFR), which all match recent observational data from the local universe. We also present high-redshift predictions of cold gas diameters and the TFR, both of which appear to evolve markedly with redshift.