We present a new three-dimensional thermomechanically coupled ice sheet model of the northern hemisphere to reconstruct the Quaternary ice sheets during the last glacial cycle. The model includes basal sliding, internal representations of the surface mass balance, glacial isostasy, and a treatment for marine calving. The time dependent forcing consists of temperature and precipitation anomalies from the UKMO GCM scaled to the GRIP ice core ∂18O record. Model parameters were chosen to best match geomorphological inferences on maximum LGM extent and global eustatic sea level change. For our standard run we find a maximum ice volume of 57 x 106 km3 at 18.5 ka cal BP. This corresponds to a eustatic sea level lowering of 110 m after correction for hydro-isostatic displacement and anomalous ice resulting from defects in the PMIP climatic forcing. Of this 110 m, 82 m was stored in the North American ice sheet and 25 m in the Eurasian ice sheet. We determine the qualitative and quantitative response of the model from a comprehensive sensitivity study in which the most important parameters were varied over their respective ranges of uncertainty. Model outputs comparable to the observational record were explored in detail as a linear function along the axes of parameter space of the reference model. The method reveals the dominance of climate uncertainty when modelling the LGM configuration of the northern hemisphere ice sheets, but also highlights the role of ice rheology and basal processes for aspect ratio, and glacial isostasy and calving for the timing of maximum ice volume.