Plastic deformation in polycrystals is governed by the interplay between intra-granular slip and grain boundary-mediated plasticity. However, while the role played by bulk dislocations is relatively well-understood, the contribution of grain boundaries (GBs) has only recently begun to be studied. GB plasticity is known to play a key role along with bulk plasticity under a wide range of conditions, such as dynamic recovery, superplasticity, severe plastic deformation , etc., and developing models capable of simultaneously capturing GB and bulk plasticity has become a topic of high relevance. In this paper we develop a thermodynamically-consistent polycrystal plasticity model capable of simulating a variety of grain boundary-mediated plastic processes in conjunction with bulk dislocation slip. The model starts from the description of a single crystal and creates lattice strain-free polycrystalline configurations by using a specially-designed multiplicative decomposition developed by the authors. This leads to the introduction of a particular class of geometrically necessary dislocations (GND) that define fundamental GB features such as misorientation and inclination. The evolution of the system is based on an energy functional that uses a non-standard function of the GND tensor to account for the grain boundary energy, as well as for the standard elastic energy. Our implementation builds on smooth descriptions of GBs inspired on diffuse-interface models of grain evolution for numerical convenience. We demonstrate the generality and potential of the methodology by simulating a wide variety of phenomena such as shear-induced GB sliding, coupled GB motion, curvature-induced grain rotation and shrinkage, and polygonization via dislocation sub-grain formation.