Grand canonical Monte Carlo (GCMC) simulations of H2 sorption were performed in Mg-MOF-74, a metal-organic framework (MOF) that displays very high H2 sorption affinity. Experimental H2 sorption isotherms and isosteric heats of adsorption (Qst) values were reproduced using a general purpose materials sorption potential that includes many-body polarization interactions. In contrast, using two models that include only charge-quadrupole interactions failed to reproduce such experimental measurements even though they are the type normally employed in such classical force field calculations. Utilizing the present explicit polarizable model in GCMC simulation resulted in a Mg2+-H2 distance of 2.60 Å, which is close to a previously reported value that was obtained using electronic structure methods and comparable to similar experimental measurements. The induced dipole distribution obtained from simulation assisted in the characterization of two previously identified sorption sites in the MOF: the Mg2+ ions and the oxido group of the linkers. The calculated two-dimensional quantum rotational levels for a H2 molecule sorbed onto the Mg2+ ion were in good agreement with experimental inelastic neutron scattering (INS) data. Although the H2-metal interaction in MOFs may be thought of as a quantum mechanical effect, this study demonstrates how the interaction between the sorbate molecules and the open-metal sites in a particular highly sorbing MOF can be captured using classical simulation techniques that involve a polarizable potential.