The osmotic pressure of salt-free polyelectrolyte solutions is studied using Monte Carlo simulations. The polymer molecules are modeled as freely jointed chains of charged hard spheres, the counterions are modeled as charged hard spheres, and the solvent is a dielectric continuum. In dilute solutions, the dominant contribution to the excess part of the osmotic pressure comes from electrostatic interactions, resulting in an osmotic coefficient that decreases with increasing concentration. In concentrated solutions, the hard sphere contribution is dominant, and the osmotic coefficient is an increasing function of concentration. By considering different terms in the pressure equation, the excess part of the osmotic coefficient is decomposed into polymer-polymer, polymer-counterion, and counterion-counterion contributions. (The electrostatic contribution can only be decomposed into two contributions because of charge neutrality.) The polymer -counterion contribution is dominant over polymer or counterion contributions in both dilute and semidilute solutions. The simulations are used to test various liquid state theories for the volumetric properties of polyelectrolyte solutions. Although all the theories considered provide a qualitatively correct picture, they tend to underestimate the electrostatic contribution, thus overestimating the osmotic pressure in dilute solutions.