Phys. Rev. B 67, 195411 (2003) We have used the Path Integral Monte Carlo method to simulate a monolayer of molecular hydrogen on graphite above 1/3 submonolayer coverage. We find that at low temperature and as the coverage increases the system undergoes a series of transformations starting from the $\sqrt{3} \times \sqrt{3}$ commensurate solid near 1/3 coverage. First, a phase is formed which is characterized by a uniaxially compressed incommensurate solid with additional mass density modulations along the same direction which can be viewed as an ordered domain-wall solid with a characteristic domain-wall distance which depends on the surface coverage. At low temperature and higher coverage there is a transition to an incommensurate solid which is rotated relative to the substrate commensurate lattice. As a function of temperature the domain-wall ordering first melts into a fluid of domain-walls and at higher temperature the solid melts into a uniform fluid. Regardless of the large amplitude of quantum fluctuations, these phase transitions are analogous to those in classical monolayer films. Our calculated values for the surface coverage and temperature where these transitions occur, the calculated structure factors and specific heat are in general agreement with the available experimental results with no adjustable parameteters.