This study developed a numerical method using a high-resolution digital elevation model (DEM) to enhance the understanding of aquifer–creek interaction in a heterogeneous aquifer system. This methodology was applied to the Stoney Creek Watershed in the province of British Columbia, Canada. The interaction between the aquifer and Stoney Creek was modeled using the 3D MODFLOW-USG code. The calibration was performed using water level measurements collected from six cross-sections along the creek from 2015 to 2017, and static water levels recorded in additional wells. In the steady-state condition, the groundwater contribution to Stoney Creek ( \({D}_{\text{g}}\) ) is estimated to be 2.47 × 107 m3/year along the creek. Under transient conditions, the simulated \({D}_{\text{g}}\) showed monthly variations, with the greatest contribution occurring in January 2015 at 2.1 × 106 m3/month. Following May 2016, the \({D}_{\text{g}}\) rapidly increased, mainly as a result of the increased precipitation occurring from May to July. Through defined management scenarios, the calibrated model was used as a predictive tool to simulate the cumulative effects of influential factors on aquifer–creek interaction. The simulation results show that the \({D}_{\text{g}}\) decreased by 21% in 2026 and by 25% in 2035 when the simulation reached steady state. Also, increasing groundwater withdrawal under dry conditions caused a potential reversal of flow from the creek into the aquifer ( \({I}_{\text{sw}}\) ), where the average \({I}_{\text{sw}}\) increased from 9 × 105 m3/year in 2017 to 1.04 × 107 m3/year by 2026. The study results and the provided spatial interaction mapping along Stoney Creek could help decision-makers better manage water resource challenges in the watershed.
