Abstract: Supraglacial ponds play a significant role in the mass loss of many debris-covered glaciers in the Himalaya. Glacier surface topography and debris flux conditions are thought to govern supraglacial ponding and drainage. Existing studies, however, have not adequately investigated the relationships and feedbacks between meltwater production, debris transport, topographic evolution and ponding, because field measurements are limited in time and space, and most existing models either neglect these processes or use oversimplified assumptions. Such limitations restrict our understanding of supraglacial hydrology and introduce uncertainties in our assessments of glacier sensitivity to climate forcing. This study develops a more comprehensive numerical model to provide insights into the couplings between topographically-controlled surface ablation, meltwater drainage, ponding, and gravitational debris transport under radiative forcing. We investigate supraglacial ponding and drainage dynamics in response to different topographic and debris flux conditions through numerical simulations based on Baltoro Glacier in the Karakoram and several hypothetical scenarios. Results suggest that: 1) Supraglacial ponds make a significant contribution to the total ice loss (more than 20 %) in the lower-mid ablation zone over one ablation season, which elevates the glacier’s nonlinear response to radiative forcing. 2) Gravitational debris transport has a non-negligible control on the growth rate of supraglacial ponds by governing debris thickness and ablation rates on the ice-cliffs around ponds. 3) Glacier surface gradient and local topographic depressions control pond formation by affecting supraglacial water storage and drainage. Our simulations provide a possible explanation to the abundance of ponds in the mid ablation zone where slope is gentle and more local depressions are present. These findings may contribute to more accurate predictions of future glacier changes in response to climate change.