Classical molecular dynamics simulations of quantum crystals, using a simple pseudopotential, are reported. The method is implemented for calculating both equilibrium and dynamical properties. As a test, the radial distribution functions for pure solid H-2 and Li doped solid H-2 are computed, and found to be in excellent agreement with prior results obtained by the variational Einstein model [D. Li and G. A. Voth, J. Chem. Phys. 96, 5340 (1992)]. The method also yields a realistic phonon density of states, which is obtained by normal mode analysis. As an implementation in dynamics, the rotational motions of O-2 isolated in solid D-2 are investigated. In agreement with recent experiments, it is found that O-2 does not rotate in solid D-2.