Simulations of the most recent 100 kyr glaciation-deglaciation cycle of the late Pleistocene Ice Age are presented which feature a newly constructed thermomechanically coupled three-dimensional ice sheet model which is itself coupled to a previously employed global energy balance climate model with realistic geography. The model incorporates both orbital insolation forcing due to the slow time evolution of orbital geometry (arising from many body effects in the solar system) as well as the forcing due to varying atmospheric concentrations of greenhouse gases. Simulations of the Greenland ice sheet are presented with which we are able to investigate the extent to which the ice flow law employed in the ice dynamics component of the model is constrained. The good agreement that we are able to achieve between the model-generated ice sheet topography and the observed Greenland topography provides a clear demonstration of the quality of the ice dynamical model we have developed. The incorporation of full thermomechanical coupling using the standard Glen flow law in the ice dynamics component of the model is shown to increase the difficulty of achieving complete termination of the 100 kyr cycle that the model delivers. Furthermore, a 20 fold flow parameter enhancement relative to that used for Greenland is required to match the aspect ratio of the ICE-4G reconstruction <Peltier, 1994>. Analyses presented herein suggest that basal processes are unlikely to account for this need for flow parameter retuning. However, the new thermomechanically coupled model now provides a clear separation of the Cordilleran and Hudson Bay domes at Last Glacial Maximum in contradistinction to analyses previously performed with isothermal ice sheet models. This innovation therefore leads to a major improvement of the generated ice sheet topography in relation to geological inferences. Aside from this important difference, the overall results regarding the ability of the model to fully account for the most recent 100 kyr cycle of glaciation and deglaciation without the necessity of introducing additional ad hoc feedbacks confirms the validity of this conclusion, reached previously on the basis of isothermal model integrations. ¿ 1999 American Geophysical Union