Aims. We aim to investigate the abundances of light deuterium-bearing species such as HD, H 2 D + , and D 2 H + in a gas-grain chemical model that includes an extensive description of deuterium and spin-state chemistry, in physical conditions appropriate to the very centers of starless cores. Methods. We combined a gas-grain chemical model with radiative transfer calculations to simulate density and temperature structure in starless cores. The chemical model includes new reaction sets for both gas phase and grain surface chemistry, including deuterated forms of species with up to four atoms and the spin states of the light species H 2 , H + 2 , and H + 3 and their deuterated forms. Results. We find that in the dense and cold environments attributed to the centers of starless cores, HD eventually depletes from the gas phase because deuterium is efficiently incorporated into grain-surface HDO, resulting in inefficient HD production on grains for advanced core ages. HD depletion has consequences not only on the abundances of, e.g., H 2 D + and D 2 H + , whose production depends on the abundance of HD, but also on the spin state abundance ratios of the various light species, when compared with the complete depletion model where heavy elements do not influence the chemistry. Conclusions. While the eventual HD depletion leads to the disappearance of light deuterium-bearing species from the gas phase on a relatively short timescale at high density, we find that at late stages of core evolution, the abundances of H 2 D + and D 2 H + increase toward the core edge, and the distributions become extended. The HD depletion timescale increases if less oxygen is initially present in the gas phase, owing to chemical interaction between the gas and the dust preceding the starless core phase. Our results are greatly affected if H 2 is allowed to tunnel on grain surfaces, and therefore more experimental data is needed not only on tunneling but also on the O + H 2 surface reaction in particular.