Mechanical twinning is a form of inelastic deformation observed extensively in low symmetry crystals, such as magnesium and other hexagonal close-packed (HCP) metals. The presence and formation of twins have been shown to have drastic effects on the strength and ductility of these metals. With recent increased interest in using HCP metals in structural, automotive, and armor applications due to their high strength to weight ratio, there is a need for comprehensive modeling of twinning to understand these material’s responses. Past studies on mechanical twinning have taken either a microscopic approach, through molecular dynamic simulations to explore sub-grain physics, or a macroscopic approach, using simplified models such as “pseudo-slip” to obtain bulk responses. However, twins interact across the mesoscale, forming collectively across multiple grains with complex local morphology propagating into bulk behavior. These mesoscale aspects have been less studied in the presence of HCP’s slip systems and so is the focus of the current work. We propose a mesoscale model where twinning is treated using a phase-field approach, while dislocation slip is considered using crystal plasticity. Lattice reorientation, length-scale effects, interactions between dislocations and twin boundaries, as well as twin and slip interactions with grain boundaries are all carefully accounted for within the model. We present the model, the implementation using a novel approach of accelerated computational micromechanics and demonstrate three-dimensional simulations of polycrystalline materials. We summarize the insights gained from these studies and the implications on the macroscale behavior of HCP materials.