Subglacial drilling through polar ice–rock transitional zones is frequently hindered by severe load fluctuations and instability induced by strong formation heterogeneity. This study investigates the coupled drilling mechanics across pure ice, an ice–rock mixture, and bedrock using a fracture-mechanics-based finite–discrete element method (F-DEM) with a strain-softening constitutive model for brittle ice. The simulations quantify the evolution of drilling force components and torque under varying gravel content, weight on bit (WOB), and rotational speed. Results show that drilling responses are strongly governed by formation-dependent mechanical behavior, transitioning from stable periodic fluctuations in pure ice to impact-dominated irregular variations in mixtures, and finally to high-level continuous resistance in bedrock. Rotational speed significantly influences drilling stability through its interaction with formation properties; lower speeds effectively suppress transient load fluctuations in pure ice and mixtures, whereas moderate speeds mitigate stick–slip effects in rock. Conversely, WOB primarily controls axial penetration while also influencing load fluctuation characteristics. The optimal WOB is identified as 10 kN for stable shearing in pure ice, reduced to 7.5 kN in mixtures to buffer lateral impacts, and increased to 10 kN in rock for efficient penetration and stable cutting. Furthermore, increasing gravel content amplifies force and torque sensitivity, transitioning the system from stable periodic behavior to impact-driven instability. These findings emphasize the necessity of a layer-dependent adaptive control strategy, providing a quantitative framework to optimize parameters and drill-bit design for safer, more efficient subglacial operations.