A model for molecular hydrogen interacting with aliphatic and aromatic hydrocarbons is presented. The model has been derived using ab initio techniques and molecular dynamics simulations. In particular, quadrupole moments of hydrogen, and variation on energy with intermolecular distance of different conformations for the hydrogen-benzene couple were calculated using the Møller-Plesset method. Hydrogen was modelled using a two-centre Lennard-Jones potential plus electrostatic interactions. Lennard-Jones parameters were optimized on the basis of a correct reproduction of experimental data of hydrogen solubility in benzene and cyclohexane, calculated using the test particle insertion method. Different sets of parameters for specific interactions (hydrogen-aliphatic and hydrogen-aromatic systems) were considered avoiding the simple use of Lorentz-Berthelot combining rules. Additionally, structural and thermodynamic properties of hydrogen-benzene, hydrogen-cyclohexane and hydrogen in an equimolar mixture of benzene-cyclohexane at different low concentrations of hydrogen were investigated by means of molecular dynamics simulations. Electrostatic charges were taken from ab initio quantum mechanical calculations but after careful analysis of the calculated properties, their irrelevance was evidenced. Moreover, Coulombic interactions make simulations more expensive and, therefore, we do not recommend their inclusion in the modelling of hydrogen-aliphatic and aromatic interactions.