Millennial-scale variability associated with Dansgaard-Oeschger events is arguably one of the most puzzling climate phenomena ever discovered in paleoclimate archives. Here, we set out to elucidate the underlying dynamics by conducting a transient global hindcast simulation with a 3-D intermediate complexity earth system model covering the period 50 to 30 ka BP. The model is forced by time-varying external boundary conditions (greenhouse gases, orbital forcing, and ice-sheet orography and albedo) and anomalous North Atlantic freshwater fluxes, which mimic the effects of changing northern hemispheric ice volume on millennial timescales. Together these forcings generate a realistic global climate trajectory, as demonstrated by an extensive model/paleo data comparison. Our results are consistent with the idea that variations in ice-sheet calving and subsequent changes of the Atlantic Meridional Overturning Circulation were the main drivers for the continuum of glacial millennial-scale variability seen in paleorecords across the globe.