Monte Carlo simulations of mesogenic polymers have been performed for single and binary systems in the bulk and for pure components in constrained environments. An expanded ensemble formalism was employed to facilitate fluctuations of the number of particles in grand canonical, Gibbs, and pseudo-Gibbs ensemble simulations. Phase diagrams were generated for pure systems of athermal 8-mers and 16-mers with varying degrees of flexibility, for a binary mixture of rigid and semiflexible 16-mers (at constant pressure and variable flexibility), and for a bidisperse mixture of rigid chains. Simulations of the partitioning of rigid 12-mer chains between a bulk and a slit pore revealed that, for a certain range of slit gaps, the preferential coplanar orientation of molecules at the proximity of the walls induces the formation of a nematic phase within the pore. A similar induction of a nematic phase was observed for anisotropic molecules dispersed in a stretched, swollen gel. For the bulk systems, our simulation results were compared to predictions of various theories,