A gradient-enhanced crystal plasticity model is presented that explicitly accounts for the evolution of the densities of geometrically necessary dislocations (GNDs) on individual slip systems of deforming crystals. The GND densities are fully coupled with the crystal slip rule. Application of the model to two distinct and technologically important crystal types, namely hcp Ti and ccp Ni, is given. For the hcp crystals, slip is permitted with a -type slip directions on basal, prismatic and pyramidal planes and c + a -type slip directions on pyramidal planes. First, a single crystal under four-point bending is simulated as the uniform strain gradient expected in the central span provides a good validation of the code. Then, uniaxial deformation of a model near- α Ti polycrystal has been analysed. The resulting distributions of GND densities that develop on the various slip system types have been compared with independent experimental observations. The model predicts that GND density on the c + a systems is approximately an order of magnitude lower than that for a -type systems in agreement with experiment. For the ccp case, slip is considered to take place on the <110>{111} slip systems. Thermal loading of a single-crystal nickel alloy sample containing carbide particles of size approximately 30 μm has been analysed. Detailed comparisons are presented between model predictions and results of high-resolution electron backscatter diffraction (EBSD) measurements of the micro-deformations, lattice rotations, curvatures and GND densities local to the nickel–carbide interface. Qualitatively, good agreement is achieved between the coupled and decoupled model elastic strains with the EBSD measurements, but lattice rotations and GND densities are quantitatively well predicted by the coupled crystal model but are less well captured by the decoupled model. The GND coupling is found to lead to reduced lattice rotations and plastic strains in the region of highest heterogeneity close to the Ni matrix/particle interface, which is in agreement with the experimental measurements. The results presented provide objective evidence of the effectiveness of gradient-enhanced crystal plasticity finite element analysis and demonstrate that GND coupling is required in order to capture strains and lattice rotations in regions of high heterogeneity.