The concept of nuclear fuel waste disposal underground is drawing increasing attention due to its many advantages against the current storage methods at surface. In this paper, we employ the Galerkin finite-element technique to solve the coupled time-dependent heat transfer and fluid flow differential equations and to predict the evolving behaviour of groundwater flow and subsurface temperature distribution associated with a proposed disposal system at the Whiteshell Research Area in southeastern Manitoba. A two-dimensional (2-D) numerical model is conceptualized from geological constraints in this particular area. To investigate the free convection of groundwater flow driven by nuclear fuel waste heat generation, we assume the 2-D model has a flat upper boundary so as to eliminate the effect of topography head. Buoyancy force due to fluid density variations is, therefore, the sole driving mechanism for fluid migration. Case studies for both unfractured and fractured porous media confirm that thermal decay of the buried fuel waste can initiate fluid circulation. In the presence of discrete fractures, the deep and hot fluid nearby the disposed contaminant can discharge to the biosphere, thus potentially threatening human health and the natural environment.