UBC Theses and Dissertations

UBC Theses Logo

UBC Theses and Dissertations

Cryohydrogeology of a covered waste rock pile in a permafrost environment : large scale field experiment… Collette, Laurier 2017

Your browser doesn't seem to have a PDF viewer, please download the PDF to view this item.

Item Metadata


24-ubc_2018_february_collette_laurier.pdf [ 20.01MB ]
JSON: 24-1.0362380.json
JSON-LD: 24-1.0362380-ld.json
RDF/XML (Pretty): 24-1.0362380-rdf.xml
RDF/JSON: 24-1.0362380-rdf.json
Turtle: 24-1.0362380-turtle.txt
N-Triples: 24-1.0362380-rdf-ntriples.txt
Original Record: 24-1.0362380-source.json
Full Text

Full Text

CRYOHYDROGEOLOGY OF A COVEREDWASTE ROCK PILE IN A PERMAFROSTENVIRONMENTLARGE SCALE FIELD EXPERIMENT ANDFREEZE-THAW NUMERICAL INVESTIGATIONSbyLaurier ColletteB.Eng., E´cole Polytechnique de Montre´al, 2015A THESIS SUBMITTED IN PARTIAL FULFILLMENT OFTHE REQUIREMENTS FOR THE DEGREE OFMASTER OF APPLIED SCIENCEinThe Faculty of Graduate and Postdoctoral Studies(Geological Engineering)THE UNIVERSITY OF BRITISH COLUMBIA(Vancouver)December 2017c© Laurier Collette 2017AbstractA large scale covered waste rock test pile (120x80x14m high) was constructedand instrumented in 2006 at the Diavik Diamond Mine in the NorthwestTerritories, Canada, in a region of continuous permafrost. The thermalcover of the pile consists of a 1.5-m thick loose till layer capped with 3 m ofnon-acid generating waste rock.The thermo-hydrologic behaviour of the core and the cover was assessedthrough the integration of meteorological, moisture content (from TDR andECH2O probes), thermistor and outflow data collected over a 10-year pe-riod. Heating cables located at the base of the pile significantly affected thethermal regime before they were turned off in 2011. The results consideredover the 10 year monitoring period show that the cover provides insulationand promotes the onset of a frozen barrier to moisture flow at the base ofthe active layer, between 2.5 and 3.5 m depth. At some locations below thecrest and within the batters, data show evidence of moisture build-up overtime on top of the frozen zone.One-dimensional freeze-thaw numerical simulations of coupled heat trans-fer and moisture flow were carried out, using SoilVision SVFlux and SVHeatpackages. The results indicate that the active layer is contained within theupper waste rock layer, where moisture builds up over time. Air temper-ature and material thermal properties are the key factors controlling thecover behaviour. The initial moisture conditions and hydraulic propertiesof the cover materials are important in the prediction of ice and moisturecontent within the till. Simulations show drainage of moisture from thetill to the waste rock in the early stage of the pile, suggesting that the tilliilayer in the test pile may not be functioning as a high moisture and latentheat layer. A long term simulation suggests that evaporation, precipitation,moisture content and the active layer thickness eventually reach a dynamicequilibrium after over 30 years if air temperature increase is neglected.iiiLay SummaryWaste rock generated by mining activities is typically disposed into stock-piles. When the rock contains sulphide minerals, rainwater infiltratingthrough stockpiles can become acidic. This phenomenon, known as acidmine drainage, is a challenging environmental issue. A common way to pre-vent acid mine drainage is to control water infiltration through a protectivesoil cover. However in regions of permafrost, the performance of covers iscomplicated by freeze-thaw dynamics.This study characterizes the hydrologic and thermal behaviour of a cov-ered waste rock pile located in a permafrost environment, based on fielddata interpretation and numerical simulations. The results indicate that afrozen layer forms within the cover and acts as a flow barrier, above whichmoisture accumulates over time to form a saturated layer. These processeswere found to be strongly controlled by thermal properties of the cover andclimatic conditions.ivPrefaceThis thesis is original and unpublished work. Unless stated otherwise, all thework including data interpretation and numerical modelling was conductedby the author, under the guidance from Leslie Smith. All chapters werewritten by the author and reviewed by Leslie Smith.The Diavik Waste Rock Project is a collaborative research program in-volving professors, graduates students, coop students and lab techniciansfrom the University of British Columbia, the University of Alberta, the Uni-versity of Waterloo and Carleton University. Most of the raw datasets thatthis thesis is based on were collected by the research team prior to the in-volvement of the author. The author’s contribution to the project includeson-going sample collection and instrument maintenance during field season2016.vTable of ContentsAbstract . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . iiLay Summary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . ivPreface . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . vTable of Contents . . . . . . . . . . . . . . . . . . . . . . . . . . . . viList of Tables . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . ixList of Figures . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . xAcknowledgements . . . . . . . . . . . . . . . . . . . . . . . . . . . xii1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 11.1 Background and significance . . . . . . . . . . . . . . . . . . 11.1.1 Cryohydrogeology . . . . . . . . . . . . . . . . . . . . 11.1.2 Numerical modelling . . . . . . . . . . . . . . . . . . 91.1.3 Cold regions . . . . . . . . . . . . . . . . . . . . . . . 91.1.4 Waste rock . . . . . . . . . . . . . . . . . . . . . . . . 101.2 The Diavik Waste Rock Project . . . . . . . . . . . . . . . . 131.3 Scope and objectives . . . . . . . . . . . . . . . . . . . . . . 161.4 Thesis organization . . . . . . . . . . . . . . . . . . . . . . . 172 Field Observations. . . . . . . . . . . . . . . . . . . . . . . . . . 182.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . 182.2 Experimental site . . . . . . . . . . . . . . . . . . . . . . . . 202.2.1 Site conditions . . . . . . . . . . . . . . . . . . . . . . 202.2.2 Test Piles . . . . . . . . . . . . . . . . . . . . . . . . . 212.3 Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 242.3.1 Thermal regime . . . . . . . . . . . . . . . . . . . . . 242.3.2 Flow regime . . . . . . . . . . . . . . . . . . . . . . . 35vi2.4 Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . . 512.4.1 Thermal barrier . . . . . . . . . . . . . . . . . . . . . 512.4.2 Barrier to flow and moisture build up . . . . . . . . . 512.4.3 Thermal influence of external heating source . . . . . 542.4.4 Heterogeneity . . . . . . . . . . . . . . . . . . . . . . 552.4.5 Proposed conceptual model . . . . . . . . . . . . . . . 562.5 Conclusion and recommendations . . . . . . . . . . . . . . . 593 Numerical Investigations. . . . . . . . . . . . . . . . . . . . . . 613.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . 613.2 Modelling . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 633.2.1 Numerical implementation . . . . . . . . . . . . . . . 633.2.2 Model Geometry . . . . . . . . . . . . . . . . . . . . . 673.2.3 Modelling assumptions . . . . . . . . . . . . . . . . . 673.2.4 Parametrization . . . . . . . . . . . . . . . . . . . . . 693.2.5 Boundary conditions . . . . . . . . . . . . . . . . . . 753.3 Base Case Results . . . . . . . . . . . . . . . . . . . . . . . . 783.3.1 Heat . . . . . . . . . . . . . . . . . . . . . . . . . . . 783.3.2 Flow and moisture . . . . . . . . . . . . . . . . . . . . 803.4 Sensitivity analysis . . . . . . . . . . . . . . . . . . . . . . . 843.4.1 Element size . . . . . . . . . . . . . . . . . . . . . . . 853.4.2 Thermal conditions . . . . . . . . . . . . . . . . . . . 863.4.3 Thermal Properties . . . . . . . . . . . . . . . . . . . 903.4.4 Flow conditions . . . . . . . . . . . . . . . . . . . . . 943.5 Long term atmospheric coupling . . . . . . . . . . . . . . . . 993.5.1 Numerical implementation . . . . . . . . . . . . . . . 993.5.2 Climatic parametrization . . . . . . . . . . . . . . . . 1003.5.3 Results . . . . . . . . . . . . . . . . . . . . . . . . . . 1013.5.4 Climate warming . . . . . . . . . . . . . . . . . . . . 1033.6 No cover simulation . . . . . . . . . . . . . . . . . . . . . . . 1043.6.1 Model set-up . . . . . . . . . . . . . . . . . . . . . . . 1043.6.2 Results . . . . . . . . . . . . . . . . . . . . . . . . . . 1053.7 Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1053.7.1 Integration of results . . . . . . . . . . . . . . . . . . 1053.7.2 Implications for cover design and construction . . . . 1103.7.3 Comments on freeze-thaw modelling . . . . . . . . . . 1123.8 Conclusion . . . . . . . . . . . . . . . . . . . . . . . . . . . . 115vii4 Conclusion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1174.1 Key findings . . . . . . . . . . . . . . . . . . . . . . . . . . . 1174.2 Recommendations . . . . . . . . . . . . . . . . . . . . . . . . 1204.3 Concluding remarks . . . . . . . . . . . . . . . . . . . . . . . 121References . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 122viiiList of Tables2.1 Covered pile materials . . . . . . . . . . . . . . . . . . . . . . 222.2 ECH2O probes thermal data summary . . . . . . . . . . . . . 352.3 Increase of the annual minimum and annual maximum VMC 423.1 Base Case hydraulic properties . . . . . . . . . . . . . . . . . 733.2 Base case volumetric heat capacity (Pham, 2013) . . . . . . . 733.3 Approximate node spacing . . . . . . . . . . . . . . . . . . . . 853.4 Sensitivity analysis results for thermal conditions . . . . . . . 893.5 Sensitivity analysis results for thermal properties . . . . . . . 913.6 Sensitivity analysis results for SFCC . . . . . . . . . . . . . . 943.7 Sensitivity analysis results for SWCC and initial moisture . . 953.8 Sensitivity analysis results for the till SWCC custom fit . . . 963.9 Sensitivity analysis results for saturated hydraulic conductivity 983.10 Climatic parametrization . . . . . . . . . . . . . . . . . . . . 100ixList of Figures1.1 Example of SFCC of different silts . . . . . . . . . . . . . . . 31.2 Hydraulic conductivity functions of partially frozen silts . . . 41.3 Temperature profile in a perennially frozen soil . . . . . . . . 51.4 Typical active layer profile in a permafrost environment . . . 61.5 Freezing and thawing indices . . . . . . . . . . . . . . . . . . 71.6 Examples of field observations of the zero-curtain effect . . . 81.7 Climate types and covers options (After INAP, 2017) . . . . . 121.8 Thermal Cover concepts . . . . . . . . . . . . . . . . . . . . . 121.9 Site location and geographical permafrost distribution . . . . 131.10 Aerial view of the experimental waste rock piles . . . . . . . . 152.1 Location map . . . . . . . . . . . . . . . . . . . . . . . . . . . 202.2 Experimental covered pile set-up . . . . . . . . . . . . . . . . 232.3 Type 3 core thermistor strings temperature profiles . . . . . . 262.4 Temperature profile of the west borehole at chainage 0+094 . 272.5 Annual temperature swing in covered and uncovered piles . . 282.6 Active layer thickness at chainage 0+094 . . . . . . . . . . . . 292.7 Close up on some Type 3 core thermistor strings . . . . . . . 302.8 Type 3 batter thermistor strings temperature profiles . . . . . 342.9 Rainfall and infiltration data . . . . . . . . . . . . . . . . . . 362.10 Outflow from the basal drain . . . . . . . . . . . . . . . . . . 372.11 C2E2 - Moisture content from TDR probes . . . . . . . . . . 382.12 C2W2 - Moisture content from TDR probes . . . . . . . . . . 392.13 C3W2 - Moisture content from TDR probes . . . . . . . . . . 402.14 ECH2O probe results at chainage 0+094 . . . . . . . . . . . . 442.15 Unfrozen moisture content, top of till layer at STA 0+094 . . 462.16 ECH2O probe results on the East batter E10 . . . . . . . . . 462.17 VMC patterns upon thawing at E10 top of till probe . . . . 482.18 VMC behaviour within the batter and under the crest . . . . 492.19 Bottom of till unfrozen moisture content at E10 . . . . . . . . 502.20 ECH2O Results in the west batter . . . . . . . . . . . . . . . 50x2.21 Integrated thermal hydrologic conceptual model . . . . . . . . 573.1 Soil water characteristic curves (SWCC) . . . . . . . . . . . . 713.2 Base case frozen hydraulic conductivity of till and waste rock 753.3 Base case model boundary conditions . . . . . . . . . . . . . 763.4 Modelled temperature compared to field data . . . . . . . . . 793.5 Base case cumulative flux (Negative flux means downward flow) 803.6 Base case results: saturation through each model region . . . 813.7 Base Case thermo-hydrology results . . . . . . . . . . . . . . 823.8 Cumulative distribution of maximum daily fluxes . . . . . . . 843.9 Element size effects on hydraulic parameters and water balance 873.10 Modified boundary conditions for temperaure . . . . . . . . . 893.11 Thermal conductivity definition - case 20 vs. base case . . . . 923.12 SFCC based on SWCC and defined as linear functions . . . . 933.13 Effects of SWCC and intial conditions on till layer drainage . 973.14 Simulated infiltration fit to data . . . . . . . . . . . . . . . . 1013.15 Atmospheric coupling results . . . . . . . . . . . . . . . . . . 1023.16 Atmospheric coupling: Active layer thickness . . . . . . . . . 1033.17 Bottom temperature condition for the cover-free model . . . . 1043.18 Dispersion of the results . . . . . . . . . . . . . . . . . . . . . 1083.19 Moisture and ice behaviour vs. temperature upon thawing . . 114xiAcknowledgementsThank you Leslie Smith for providing judicious advices, for your sense ofhumour and for the confidence towards the advancement of my work. Mod-elling in a research context can go in many directions and face many un-knowns and obstacles. I am thankful that you have always provided preciousguidance at the right time.Thanks to Jordan Zak, for sharing your knowledge about the test piles, bothin the office and at the site.Thanks to my other office friends: Jarod, Keelin, Andrea, Laura and Bren-dan. It was definitely better when you guys were all around! Thanks toConstanza and Jessie, I truly appreciated the discussions we had together.I wish you the best for the completion of your degree, and the most excitingfuture!Thanks to David Wilson and Jeff Bain for helping me to solve some of themysteries within the datasets. Thanks to all the test pile crew and theDDMI environment team.Thanks to Jirka Sˇimu˚nek for all the very quick and clear responses on theHydrus-1D forum and for recompiling the freezing module after that incon-sistencies were noted. It was much appreciated.Thanks to Nam Pham and the SoilVision technical support.Thanks to all the awesome people I met in Vancouver and BC in the lasttwo years, through all these amazing adventures... it was a blast!I finally want to thank my parents and family for the constant love, encour-agement and support... Merci infiniment pour tout!xiiFunding of this work was provided by NSERC (Discovery Grant) held byLeslie Smith and through NSERC (CGS M), Colin D Spence Memorial andCSAP (Contaminated Sites Approved Professionals Society) scholarships.Field data used in this thesis was obtained during the course of the DiavikWaste Rock Project, jointed funded by an NSERC CRD and Diavik Dia-mond Mine.xiii1Introduction1.1 Background and significance1.1.1 CryohydrogeologyIn cold regions and permafrost environments, the movement of moisture andthe occurrence of ground ice in variably saturated soils (natural or man-made) are closely related to soil freezing and thawing dynamics, influencedby thermal conditions and soil properties. The study of these interdepen-dent phenomena, refered to as cryohydrogeology, has so far received littleattention in the mining environmental literature.Key conceptsThis first section stands as a short literature review by providing brief defini-tions of the key concepts relative to permafrost environment and soil freezingand thawing. Terms relevant to this thesis are presented.Latent Heat (of fusion or of vaporization) corresponds to the amount ofenergy absorbed or released during phase change. The amount of energyabsorbed or released upon increase or decrease of temperature is rather1called sensible heat. In a given soil, the amount of latent heat involvedupon phase change is a function of the total amount of moisture in the porespace as well as the fraction of the water that experiences phase change(Andersland & Ladanyi, 1994). This fraction is linked to the Soil-FreezingCharacteristic Curve (defined below). Latent heat absorption and release inthe subsurface can be sufficient to slow down or stop the progression of thefreezing or thawing front (McFadden & Bennett, 1991).Clausius-Clapeyron Equation (CCE) expresses the relationship be-tween the pressure of ice and of liquid water in freezing soil and temperature.It is derived from the thermodynamic Gibbs free energy concept (Kurylyk& Watanabe, 2013). The original form of the CCE is the following:dPwρw− dPiρi= LfdTT(1.1)where Pw is the equilibrium gauge pressure of the liquid water phase, ρw thedensity of water, Pi is the equilibrium gauge pressure of the ice phase, ρithe density of ice, Lf the latent heat of fusion and T the temperature. Theequation is only valid at thermodynamic equilibrium. At disequilibrium,the application of the CCE tends to underestimate the unfrozen moisturecontent and hydraulic conductivities in frozen soil (Watanabe & Osada,2016).Soil-Freezing Characteristic Curve (SFCC) defines the unfrozen vol-umetric moisture content (UVMC) in a soil at a given sub-zero temperature.This curve is very similar to the Soil-Water Characteristic Curve (SWCC),2also known as retention curve. Recent theory development and research ex-perimentations support the assumption that the amount of unfrozen mois-ture content at sub-zero temperature depends on the soil properties and isindependent of the total water content (Kurylyk & Watanabe, 2013). Fig-ure 1.1 presents a SFCC example. Many studies report the similarities be-tween the SWCC and the SFCC (Bittelli, Flury, & Campbell, 2003; Black &Tice, 1989; Koopmans & Miller, 1966; Spaans & Baker, 1996; P. J. Williams,1964). A common way to define the SFCC is through empirical relationshipsthat relates the two curves or simply via the CCE to link temperatures andsuction values.Figure 1.1: Example of SFCC of different silts (Adapted from Azmatch etal., 2012)Hydraulic Conductivity of frozen soil At sub-zero temperatures, theunfrozen moisture content decreases significantly according to the SFCC andthe CCE. As a result, suction rises and the hydraulic conductivity decreasesthe same way it does in unsaturated soils (Figure 3.2). Methods for estimat-3Figure 1.2: Hydraulic conductivity functions of partially frozen silts(Adapted from Azmatch et al., 2012)ing hydraulic conductivity functions of a partially frozen soil are discussedin Azmatch et al. (2012) and Kurylyk and Watanabe (2013).Cryosuction is the movement of water under a pressure gradient causedby freezing and usually occurs at the frost front. This phenomenon is ex-plained by the CCE. When the gauge pressure in the ice phase is neglected,Equation 1.1 indicates that change in temperature induces a change in suc-tion. Behind the freezing front, the colder material has a lower unfrozen wa-ter content and thus a lower matric potential than the warmer layer aheadof the front. If the suction gradient that develops is great enough, water willtend to flow from the unfrozen material to the frozen zone at the freezingfront.Annual Temperature Swing (ATS) or the annual ground temperaturevariations or annual temperature amplitude is the difference between theannual maximum and minimum temperatures at a specific depth. As shown4on Figure 1.3, the ATS is larger near the surface and decreases at depth.It becomes negligible at what is called the depth of zero annual amplitude(Andersland & Ladanyi, 1994).Figure 1.3: Temperature profile in a perennially frozen soil (after Ander-sland & Ladanyi, 1994)Active Layer (AL) or seasonally frozen ground or annually thawed layercorresponds to the top layer of the ground where the temperature fluctu-ates above and below 0◦C throughout the year, as indicated on Figure 1.3.The terminology related to the AL dynamics is synthesized on Figure 1.4.Upward and Downward Freezing (or up-freezing and down-freezing) refer tothe direction of the freezing front propagation during the AL freeze-back(Luo, Jin, Lu¨, & Wu, 2014; Osterkamp & Romanovsky, 1997; Woo, 2012).5Figure 1.4: Typical active layer profile in a permafrost environmentUpward freezing usually accounts for less than 50% of the freeze-back whentwo-side freezing occurs (Luo et al., 2014; Osterkamp & Romanovsky, 1997).Surface N-factors express the relation between air temperature andground surface temperature. They are controlled by many site specific con-ditions such as solar radiation, wind, the occurrence of a snow cover andvegetation. When temperature is above zero, the thawing n-factor nt is de-fined as the ratio of the ground surface thawing index, Ist, to the air thawingindex, Iat:nt =IstIat(1.2)The same definition applies for the freezing n-factor nf when temperatureis below zero, but with ground surface and air freezing indices (Isf and Iaf ):nf =IsfIaf(1.3)6Freezing and thawing indices refer to the freezing and thawing degree-days,as illustrated on Figure 1.5.Figure 1.5: Example of freezing index for ground surface temperature andthawing index for air temperatureZero-Curtain Effect is a term that “refers to the effect of latent heatin maintaining temperatures near 0◦C over extended periods in freezing orthawing soils” (Outcalt, Nelson, & Hinkel, 1990). High ice or moisture con-tent prolongs the zero curtain effect and tends to limit the propagation ofthaw within the active layer because more latent heat is involved duringphase change (Carey & Woo, 1998; Outcalt et al., 1990; Shiklomanov &Nelson, 2013; Woo, 2012). Shiklomanov and Nelson (2013) noted that ob-servations of the zero-curtain effect during annual thawing are less frequentthan during freezing, as emphasized by the examples on Figure 1.6.7(a) Adapted from de Grandpré et al. (2012) (b) Adapted from Shiklomanov and Nelson(2013)(c) Adapted from Zona et al. (2016)(d) Adapted from Weismüller et al. (2011)(e) Adapted from Urban and Clow (2017)Figure 1.6: Various examples of field observations of the zero-curtain ef-fect, highlighted in light blue.81.1.2 Numerical modellingThe development of numerical models for simulating freeze-thaw, moistureflow and heat transfer in variably saturated soils (cryohydrogeology models)is an active field of research. Many codes based on different assumptions havebeen developed over the last decades (Kurylyk & Watanabe, 2013). Theyall come with their own set of limitations. For instance, Zhao et al. (2016)distinguish two classes of freezing-thawing models based on assumptions onice pressure. Hydrodynamic models that assume constant ice pressure failat predicting frost heave but are more accurate for unsaturated conditions.Conversely, frost heave models account for variable ice pressure. They wellpredict frost heave but are less accurate for unsaturated conditions.Hence, it is necessary to use a simulator that is adapted to the problembeing addressed. It is also essential to test the modelling approaches in lessdocumented cases such as flow through waste rock in a permafrost environ-ment, especially if experimental measurements are available for validation.1.1.3 Cold regionsThe understanding of freeze-thaw processes in variably saturated natural ornon-natural soils has numerous applications such as agriculture, land andecosystem management in cold regions (He et al., 2015; Iwakun, Biggar, &Sego, 2010; Luo et al., 2014; Shiklomanov & Nelson, 2013; Zhao et al., 2016)as well as artificial ground freezing involved into underground work (Anders-land & Ladanyi, 1994; Jones, 1981) and frozen soil barriers to contaminatedgroundwater (Andersland, Wiggert, & Davies, 1996; Wagner, 2013). It also9widely applies to permafrost and cold region engineering, which aims tounderstand the vulnerability of roads, railways, bridges, tunnels, buildingsand water management facilities at high latitudes (Al-Houri, Barber, Yonge,Ullman, & Beutel, 2009; Andersland & Ladanyi, 1994; Gholamzadehabol-fazl, 2015; Tan, Chen, Wu, & Yang, 2013; Wang, Wang, & Yu, 2015; Zhang,Zhang, Zhang, Chen, & You, 2016). Reliable heat transfer and water flowmodels for freezing and thawing soils are needed to assess the impact ofclimate change on permafrost degradation, active layer expansion and hy-drogeological systems (Kurylyk, MacQuarrie, & McKenzie, 2014; Walvoord& Kurylyk, 2016). Finally, mine closure plans in cold regions also includefreeze-thaw and permafrost to the design of tailing management facilities(Coulombe, Bussie`re, Coˆte´, & Garneau, 2012; Elberling, 2004; Kyhn & El-berling, 2001) and waste rock piles (Pham, 2013; Smith et al., 2013; Stevens,2016). The disposal of waste rock in a permafrost environment is the issueaddressed in this research.1.1.4 Waste rockAcid Mine DrainageMining activities, especially open pit operations, extract significant amountsof non-valuable material out of the ground. In Canada alone, the mostrecent available statistics reveal that the mining industry (metal and non-metal) produced 256 millions of tonnes of waste rock in 2008, according toNatural Resources Canada (2008). The waste rock is generally disposed inlarge stockpiles standing above ground, thus exposed to atmospheric con-10ditions. Because of the presence of water and oxygen, oxidation reactionspotentially occur in waste rock containing sulphide-bearing minerals. Theprocess forms Acid Mine Drainage (AMD), a topic well documented in theliterature (Akcil & Koldas, 2006; Amos et al., 2015; Nordstrom & Alpers,1999). Without mitigation, AMD is known to discharge low pH leachatewith dissolved metals and other toxic substances into the environment (Ak-cil & Koldas, 2006), which can result in serious ecological, human health andeconomical consequences (Azapagic, 2004; Jacobs, Lehr, & Testa, 2014).The significance of AMD loadings in the environment depends on thephysical and geochemical properties of waste rock, as well as processes thatcontrol moisture and oxygen flow through waste rock piles (Akcil & Koldas,2006; Jacobs et al., 2014; Nordstrom & Alpers, 1999). Studying and pre-dicting the hydrology of waste rock piles is important but rather complexbecause of internal structures and the heterogeneity of the material involved(Amos et al., 2015).CoversIn order to reduce water and oxygen infiltration contributing to AMD, mineclosure plans often include the construction of a cover on top of waste rockpiles. Soil cover design must be site specific (INAP, 2017) and accountfor climatic conditions as shown on Figure 1.7. Thermal covers are in apermafrost environment. As illustrated on Figure 1.8, a number of thermalcover concepts exists, depending on material availability and feasibility.11Figure 1.7: Climate types and covers options (After INAP, 2017)Figure 1.8: Thermal Cover concepts (Adapted from Pham, 2013 andStevens, 2016)121.2 The Diavik Waste Rock ProjectThe Diavik Waste Rock Project is a multi-disciplinary research programthat examines the long term hydrological, thermal, geochemical and gas-transport behaviour of waste rock piles in a continuous permafrost envi-ronment. The Diavik Diamond mine is located about 230 km south of theArctic Circle in the Northwest Territories (Figure 1.9). Based on data from2007 to 2015, the average on-site temperature is -8.8◦C. The minimum andmaximum monthly average is respectively -27.5◦C in January and 13.1◦Cin July. The annual precipitation varies between 200 and 250 mm (GolderAssociates, 1997).Figure 1.9: Site location and geographical permafrost distribution(Adapted form Stevens, 2016)13On site, three large scale waste rock experimental piles were built in 2006(Figure 1.10). The materials used to build the piles include:− Type 1 Waste Rock: Non-acid generating, Sulfur content less than0.04 wt.%, mostly granite with lesser amount of granitic pegmatiteand diabase− Type 3 Waste Rock: Potentially acid generating, Sulfur contentgreater than 0.08 wt.%, Granite with more biotite schist− Loose till: Surficial material (local soil overburden), poorly gradedsand with silt and gravelTwo test piles are uncovered (the Type 1 and Type 3 test piles). The thirdone is a Type 3 core covered pile, approximately 120 m × 80 m × 14 m high.The thermal cover design is intended to work as an air convection cover witha high moisture content layer (Figure 1.8). It consists of a 1.5-m thick tilllayer capped with 3 m of Type 1 material. All piles were instrumented withthermistors, moisture probes, gas sampling ports, moisture sampling ports,basal collection lysimeters, tensiometers, basal drainage collection systems,air permeability probes, microbiology access ports and thermal conductivityports. Further details on the experimental site design, construction andinstrumentation is provided by Smith et al. (2013), Neuner (2009), Fretz(2013) and Pham et al. (2013). The work presented in this thesis focuses onthe covered test pile.14Figure 1.10: Aerial view of the experimental waste rock piles (Adaptedfrom Fretz, 2013)151.3 Scope and objectivesThe goal of this research is to characterize the thermal and hydrological be-haviour of an experimental waste rock pile and its thermal cover composedof till and non acid generating waste rock in a permafrost environment. Thework includes the analysis and interpretation of field data collected over a10 year period as well as numerical investigations of a one-dimensional rep-resentation of the covered test pile. Following the extensive thermal charac-terization and modelling performed by Pham (2013), this research focuseson both moisture flow and freeze-thaw processes to provide insight to thethermo-hydro-geochemical processes within the test pile.The specific objectives of the research are the following:1. To assess the cover performance regarding:− the evolution of the active layer,− the conditions leading to the onset of a barrier to fluid flow,− moisture and ice distribution and accumulation within the cover,− and the differences between processes occurring at the top of thepile and in the batters.2. To document the history of flow and the thermal evolution within thecovered test pile.3. To determine and document the material properties or site conditionsthat have significant control on the cover behaviour.164. To investigate the onset of a long term dynamic equilibrium in thecover with respect to infiltration, evaporation and thermal conditions.5. To assess the potential feedback effects of moisture and flow patternsin the cover and the interior of the pile on the thermal regime.1.4 Thesis organizationThis thesis is composed of four chapters. The introductory Chapter 1 pro-vides background information relative to cryohydrogeology and details aboutthe context of this research. Chapters 2 and 3 were written as stand-alonearticles intended for publication. Chapter 2 presents a characterization ofthermal and hydrologic behaviour of the covered test pile based on datacollected for a 10 year period. Chapter 3 provides complementary work toChapter 2 in order to better understand the hydrology of the pile and thecover with respect to the internal thermal regime. Numerical simulationsof coupled heat transfer and moisture flow are detailed and interpreted. Fi-nally, Chapter 4 summarizes the key findings of this research and providescomments and recommendations for future work.172Field Observations of the Long TermThermo-Hydrological Behaviour of aCovered Waste Rock Test Pile in aPermafrost Environment2.1 IntroductionAcid mine drainage (AMD) in waste rock piles containing sulphide mineralsis well documented in the literature. Without mitigation, AMD is known todischarge low pH leachate with dissolved metals and other toxic substancesinto the environment (Akcil & Koldas, 2006), which can result in seriousecological, human health and economical consequences (Azapagic, 2004; Ja-cobs et al., 2014). The magnitude of AMD loadings to the environmentdepends on the physical and geochemical properties of waste rock, as wellas processes that control flow through waste rock and oxygen supply (Ak-cil & Koldas, 2006; Jacobs et al., 2014; Nordstrom & Alpers, 1999). Oneof the most common and accepted ways to reduce water and oxygen infil-tration through waste rock piles is the construction of a cover. Soil cover18design depends mostly on climatic conditions (INAP, 2017) and materialavailability.The objective of this chapter is to assess the hydrologic evolution of anexperimental waste rock pile with a multi-layer cover composed of till andrun-of-mine waste rock in a cold climate environment. In such an environ-ment, near and sub-zero temperatures are common. Water flow experiencesfreeze-thaw cycles and permafrost comes into play. These factors add com-plexity to understanding cover performance. The integration of meteoro-logical data, moisture content, temperature, infiltration and outflow ratesprovides a key to understanding the moisture distribution (water & ice)within the pile with respect to thermal regime and cover design.Long-term performance of soil covers in humid and temperate settingshas been well studied in the literature (Weeks & Wilson, 2005; D. J. Williams,Stolberg, & Currey, 2006). This research studies a large scale cover exper-iment in a dry and cold climate, monitored over a 10-year period. It ispart of the Diavik Waste Rock Project, a multidisciplinary research teaminvestigating the long term behaviour of waste rock piles in a cold climateenvironment based on a multiple scale field experiment. Previous work ex-amined geochemistry, hydrogeology, microbiology, gas transport and heattransfer (Amos, Blowes, Smith, & Sego, 2009; Bailey, Blowes, Smith, &Sego, 2015; Bailey et al., 2013; Chi et al., 2013; Neuner et al., 2013; Phamet al., 2013; Smith et al., 2013).192.2 Experimental site2.2.1 Site conditionsThe experimental waste rock test piles are located at the Diavik DiamondMine in the Northwest Territories, located about 230 km south of the ArcticCircle (Figure 2.1). The mine is on a 20 km2 island (Chi, 2010; Smith etal., 2013) in Lac de Gras. Under a flat tundra topography lays continuouspermafrost (90%-100%) of a few hundreds meter thick (Natural ResourcesCanada, 1995).10kmDIAVIKLAC DE GRASYellowknifeNORTHWEST  TERRITORIESFigure 2.1: Location mapThe mean annual rainfall for the 2007-2015 period was 99 mm. A re-port by Golder Associates (1997) and more recent data from EnvironmentCanada indicate that rainfall accounts for 40 to 55 % of total annual pre-20cipitation. For the same period, the average onsite temperature was -8.8◦C,with a minimum and maximum monthly average of respectively -27.5◦C inJanuary and 13.1◦C in July.2.2.2 Test PilesThe test piles site includes three large scale experimental waste rock piles.Two of them are uncovered piles while the third one is a covered pile. Thecovered test pile was constructed in 2006 and 2007, with the initial outflowoccuring in September 2007. Full description of design, construction andinstrumentation is detailed by Smith et al. (2013). This study focuses onthe covered pile. The main materials used to build the covered pile aredetailed in Table 2.1. As illustrated in Figure 2.2, the Type 3 waste rockcore (sulphur content greater than 0.08 wt. %) is approximately 10-m thick.Before construction of the cover, batters were re-sloped to a 3H:1V angle.The cover design consists of a 1.5-m thick till layer capped with 3 m ofType 1 waste rock (sulphur content less than 0.04 wt. %). Given how covermaterials were placed by heavy equipment (Smith, 2008), some thicknessvariability is expected: up to 0.5 m per layer according to as-built drawingsand field notes. The top of the pile (crest) is approximately 26 m wide at afinal elevation of 453.5 meters above sea level. The 14-m high covered pilehas a footprint approximately 120 m long by 80 m wide.The location of the instruments relevant to this study is shown onFigure 2.2. For moisture content measurement, 10 TDR probes (Time Do-main Reflectometry) are located in the Type 3 core, along tip faces 2 and 3.For moisture content and temperature, a total of 11 ECH2O probes below,21Table 2.1: Covered pile materialsType 1 WasteRockNon-acidgenerating(<0.04 wt.% S)Primarly granite withlesser amount of graniticpegmatite and diabaseParticle size from silt toboulders, dominated bysand and gravel.Bulk porosity around0.24Type 3 WasteRockPotentially acidgenerating(>0.08 wt.% S)Granite with more biotiteschistTillSurficial material (soil overburden)excavated during early mining pit developmentPoorly graded sand with silt and gravel (SP-SM)(0% clay, 11% silt, 57% sand and 33% coarser).above and within the till layer are located in the batters and the under thecrest of the pile. Temperature is also monitored using 21 thermistor strings(12 thermistors per string) located on tip faces 1 to 4, under the bottomliner, in the borehole at chainage 0+094, as well as in bedrock boreholesbelow the pile. Smith et al. (2013), Neuner (2009) and Pham et al. (2013)provide further details on the instruments, their installation and calibration.22Figure 2.2: Experimental covered pile set-up232.3 Results2.3.1 Thermal regimeThis section first examines the effect heat cables within the pile, then datafrom a borehole string near the northern limit of the pile, data from Type3 core, the response recorded in the batters and then the till layer.Heat traceA heat trace system was installed in the basal collection lysimeters (BCL)and along the basal drains during construction of the base of the pile, beforeplacement of waste rock. The heat trace was meant to ensure that moistureentering the BCL’s could drain to the instrumentation trailer in order toget a good set of water chemistry data (Smith et al., 2013). The heat tracewas turned off in mid-2011. It was also shut off temporarily between lateJuly and late August 2008. When turned on, the power control was suchthat it generated a significant amount of heat around each BCL. Accordingto the data (thermistor strings at the heat trace and in BCL waste rock),temperature were up to 20-30◦C in the BCL.As a result, thermal data reflect two major thermal periods for the first 10years of data. The first period consists of a warm phase, from the beginningof the experiment until around mid-2011 to early 2012. Then, the heattrace was turned off in June 2011 (Krentz, 2014). The second period is thefollowing cooling phase from 2012 onwards. The heat trace effect is greatestclosest to the BCL and dampens as distance increases away from the BCL.Thus, Type 3 material below the cover on tip-face 1 and 2 show the largest24effect (Figure 2.3) whereas noticeable but minor impact is observed in the0+094 boreholes (Figure 2.4).Borehole string at 0+094Annual Temperature Swing As revealed by thermistor strings at 0+094,with a cover, the annual temperature swing (ATS) or annual temperatureamplitude (Andersland & Ladanyi, 1994) at only 1.5 m depth is reduced byapproximately 50% relative to the ground surface variation (Figure 2.5a).Below the till layer (beyond 4.5 m depth) the ATS is less than 15% of theground surface ATS and it drops as low as 1-2% at the base of the pile. Incontrast, ATS at depth in the uncovered pile is never less than 40% of ATSnear the ground surface of the pile (Figure 2.5b).Changes of ATS from one year to another in the covered pile is dampenedwith depth. This is less obvious in the uncovered pile, where ATS varies alot more at any specific depth. In the uncovered pile, there is no evidenceof a relationship between the depth of measurement and the rate of changein ATS from peaks to troughs. The ATS evolution with time is more stableeverywhere in the covered pile, whether it is below, within or above the tilllayer. This indicates that there is an influence of the thermal regime in thetill layer on the waste rock forming the upper cover.Active Layer The active layer thickness varies by 12 cm on average fromone year to another (Figure 2.6), following the aforementioned warming andcooling periods. On both the east and west sides, thickness increases alongwith internal warming, from the beginning of data collection until 2012.25Figure 2.3: Type 3 core thermistor strings temperature profiles.Naming convention: C stands for covered pile, 1-4 is the face number (seeFigure 2.2), W/E stands for west or east from the centerline, 5 is the offsetfrom the centerline (in meters).To illustrate thermistors measurement uncertainty of 0.2◦C, the -0.2, 0 and0.2 isotherms are highlighted in black.26Figure 2.4: Temperature profile of the west borehole at chainage 0+094(refer to Figure 2.2). Black lines are the -0.2, 0 and 0.2 isotherms (0.2◦C isthe thermistor reading error)Later on, as the pile cools down from, the active layer thickness is reducedand stabilizes around 2.5-3 m, within the Type 1 waste rock above the tilllayer. For comparison, the active layer thickness of the Type 1 uncoveredpile was about 12 m (Zak, 2017).The ground surface temperature is above zero for a 148-day period onaverage (stdv 9 days). As there is no trend over the years, this relativelyconstant surface condition support the assumption of steady climatic con-ditions discussed later. Changes in active layer thickness as well as freezingand thawing front velocities are then most likely independent of air temper-ature fluctuations. They would be due to the internal heat source, alongwith rising of permafrost within the pile.Up-freezing (upward freezing from the bottom of the active layer) anddown-freezing (downward freezing from the ground surface) are common inlayers overlying permafrost (Luo et al., 2014; Woo, 2012). Data at 0+094show the occurrence of such bidirectionnal freezing. The average ratio ofup-freezing thickness over the total active layer thickness is 18% (varying27(a) Annual temperature swing in the covered pile (0+094 east borehole)(b) Annual temperature swing in the uncovered pile (T31S5)Figure 2.5: Annual temperature swing in (a) covered pile and (b) un-covered pile, based on monthly averaged temperature at each thermistorbead.28Figure 2.6: Active layer thickness at chainage 0+094from 7% to 29%). Therefore, for the study period, down-freezing countsfor over 2/3 of active layer freeze-back. The same ratio was observed byOsterkamp and Romanovsky (1997) in a study of a natural soil profile.Type 3 core thermistor stringsStarting from July 2010, raw data from half of the thermistors on stringslocated in the west side of the Type 3 core shifted up by few degrees. Ther-mistors located on the east side do not show this effect. The erroneous timeseries were corrected based on the last valid data points before the shiftoccurs. Since the cause of the problem remains unknown, the uncertaintycannot be assessed precisely. Corrected absolute temperatures and resultingisotherms are likely off by one or two degrees, but relative variations suchas cooling, warming and flat stable periods were found to be adequate forthe following analysis.As mentioned earlier, the heat trace in the BCL’s had a significant im-pact on the thermal regime in the Type 3 core of the pile. The strings located29Figure 2.7: Close up on some Type 3 core thermistor strings temperaturetime series. Dashed lines represent thermistors measurement uncertainty of0.2◦C. Naming convention: C stands for covered pile, 1-4 is the face number(see Figure 2.2), W/E stands for west or east from the centerline, 5 is theoffset from the centerline (in meters).30at tip face 4 show the coldest temperature and strings at face 1 and 2 showwarmer temperature. Warming and cooling periods are obvious at all loca-tions (Figure 2.3). The analysis in this chapter focuses on the cooling period(2012 onwards). With the removal of the artificial heat source, the systemis getting closer to ‘natural’ conditions. At all string locations, freezing isoccurring from the top downwards, following an annual cycle controlled bythe atmospheric temperature variations at the the pile surface.At some locations the temperature response shows a pause of the zeroisotherm that last several months (Figure 2.7). The phenomenon is referredto as the zero curtain effect in the literature (Outcalt et al., 1990; Shiklo-manov & Nelson, 2013; Sumgin, 1927; Woo, 2012). It is a result of latentheat transfer involved in phase change of water. In the present case, thiseffect is indicating that a significant amount of pore ice is slowly melting.The phenomenon is not seen at all depths on some strings, but it is at leastevident at the topmost thermistor(s). Data show that where the zero cur-tain effect occurs, it lasts for over 50% of the thawing period. These zonesof the core are not experiencing dry conditions. For instance on the westside, the topmost thermistors in the core, located 0.5 to 1 m below the till,show that temperature always remains below zero year-round (Figure 2.3).The shape of these thermistor time series (Figure 2.7 - C2W5 at 5.2 m) in-dicates frozen conditions and suggest high pore ice and water content sincethe zero-curtain effect lasts for almost half a year.The presence of ice is such that it significantly affects the evolution ofthe thermal regime. Phase change is known to affect cooling and warmingrates (Andersland & Ladanyi, 1994). In non-permafrost terrain, latent heat31effects tend to hinder quick thermal changes and thus prevent frost frontpenetration (McFadden & Bennett, 1991). Conversely, it is assumed thatwith negligible ground ice content, temperature would have risen above zeroat depths where zero curtain effects are manifest.Despite the variability between one string to the other, the general pic-ture is the same at all string locations. In thermistors at C3E5, the annualamplitude is less than 1◦C whereas it varies from 1 to 4◦C at C4E5 (Fig-ure 2.7). At C2E5, temperatures at depth during the cooling period areamong the coldest. There is an consistent frozen zone around 10 m in thatzone (Figure 2.3). Otherwise, the temperature history in the core is similarbetween the east and west sides.Type 3 batter thermistor stringsThermal data from the batter strings also show the effect of heat tracewarming and shut off. On average, temperature are warmer towards face 1and cooler towards face 4.At all of the 4 tip faces, data indicates that the waste rock below the tillon the east side is colder than on the west side. During the cooling period,the freezing front advances below the till layer on east side at all faces. Thatis not the case everywhere on west side. Indeed, unfrozen patches seen inFigure 2.8 are smaller in volume and shorter in duration on the east side.The till layer is a slightly thicker on the west side according to the as builtdrawings. This may explain the East-West difference. It is however mostlikely due to differential wind-induced cooling of the pile. Amos et al. (2009)32found a main wind direction predominance from the East and a secondarypredominance from the North.Interestingly for each string, data show warmer temperature in the mid-dle of the batter and consistently colder patterns at the toe and near thecrest. This is consistent with the observation that the crest is more affectedby the wind and the toe is closer to permafrost table. In the region betweenthe crest and the toe, snow provides insulation to the batters. The snowpack in the winter was negligible at the crest and increases towards the baseof the pile (Krentz, 2014).On the east side, during the cooling period, the farthest thermistordownslope of all the ones that experience annual thawing shows an increas-ing zero curtain effect year after year. Assuming similar annual thermalconditions (freezing and thawing rates), it suggests some moisture buildingup. Somehow, moisture cannot drain on the slope at these locations becauseflow is likely hindered by a frozen barrier.Till batter thermistor stringOnly one data set from the batter strings is available from thermistors lo-cated within the till layer. As the thermistors are closer to the surface,they show colder temperatures in the winter than the corresponding stringsbelow, in the Type 3 waste rock core. For the whole study period, all tillthermistors except the closest to the toe experience annual freeze-thaw. Thisindicates that they are within the batter active zone.Four locations (-30.8, -28.9, -27.0 -25.1 meters off the centreline) showa clear increase of the zero curtain effect in the last 4 years. More data33Figure 2.8: Type 3 batter thermistor strings temperature profiles.To illustrate thermistors measurement uncertainty of 0.2◦C, the -0.2, 0 and0.2 isotherms are highlighted in black. On each plot, two solid vertical linesshow the horizontal location of the edge of the crest, approximately 13 moff the centerline.34Table 2.2: ECH2O probes thermal data summary, stressing on tempera-ture behavior upon thawingApproximatedepthBelow zeroyear-roundSlightly abovezero (1-2◦C)Up to 5◦Cabove zeroUp to 7 to 20◦Cabove zero, withsurface temp.oscillationType 1 WR 0+094Top till 0+094 & E16 E10 W10 & W16Bottom till 0+094 & E16 E10 W10 & W16 SouthType 3 WR Southwould help to test the hypothesis that the total volumetric moisture content(VMC) is increasing in that zone.ECH2O probes temperaturesTemperature patterns previously discussed based on thermistor beads areconsistent with the ECH2O thermal data. The summary presented in Table2.2 confirms the trend that the cover on the east side is overall colder thanon the west side. It also appears that the south batter is even warmer thanthe west batter.2.3.2 Flow regimeThis section presents the hydrology of the covered pile by first introducingthe infiltration and outflow data. Then, moisture content trends in the coreare presented, and then moisture data in the cover.35InfiltrationAs described in Zak (2017), a modified Penman-Monteith (PM) model wasused to calculate the infiltration at the surface of the experimental wasterock piles. The estimates for 2008 to 2011 fit within 91-108% of in situ waterbalances (Fretz, 2013). For the 10-year study period, the annual infiltrationvaries from 9 to 85 mm (Figure 2.9). The average is 38 mm, which is 38%of the mean annual rainfall. The infiltration averaged over the whole pilelikely differ from this number because of the dominance of batter surfacewith respect to the flat crest surface. Indeed, less infiltration is expected onbatters because of run-off on the slopes.Figure 2.9: Rainfall and infiltration dataOutflowThe pile outflow follows an increasing trend from 2007 to 2011 (Figure 2.10).The following year, total annual outflow drops significantly. Indeed, theoutflow of 2012 is 7% of that of 2011. In 2013, it’s 47% of 2012 outflow. From2014 onwards, the outflow becomes negligible. Consistent with thermalregime, flow ceases as cooling initiates.36Figure 2.10: Outflow from the basal drainMoisture content in the coreAnnual wet-up and drainage The first pattern observed in the moisturecontent data from the TDR probes is a major annual wet-up (rapid VMCincrease) and drainage (VMC decrease) sequence. This is clearly demon-strated at C2E2 and C2W2 (Figures 2.11-2.12). Under unfrozen conditions,a quick wet-up occurs upon the arrival of a wetting front moving down-wards, likely triggered higher up in the pile by the thawing front reachingthe bottom of the active layer.The wet-up occurs in the shallowest probe first approximately half-waybetween the start and end of the rainfall season, 2-3 months after the firstrain event. Then, it advances from top down to the bottom of the pile.Upon wet-up, VMC reaches its maximum value for the year.It takes a relatively short time for the soil to change from a low unfrozenVMC to wet conditions (high unfrozen VMC). It occurs over roughly 10 to12 days on average. At greater depths, it can take a few weeks up to a37Figure 2.11: C2E2 - Moisture content from TDR probes on tip face 2of the covered pile (C2), 2 m east of the centerline (E2). Dashed linesemphasize multi-year trends.month. For the rest of the time, VMC is decreasing until the arrival of thenext wetting front the following year. The drainage is a much slower andmore irregular process compared to the annual wet-up.According to the portion of data showing a clear downward progressionof the wetting front, the average velocity of the wetting front is 7 cm/d(8 · 10−7 m/s) in the waste rock core below the till, with a standard devi-ation of 4 cm/d. Assuming that water is moving under a unit gradient inunsaturated material, the values are consistent with estimations of saturated38Figure 2.12: C2W2 - Moisture content from TDR probes on tip face 2of the covered pile (C2), 2 m west of the centerline (W2). Dashed linesemphasize multi-year trends.hydraulic conductivity by Neuner (2009) when accounting for moisture con-tent. Neuner (2009) found a geometric mean of 9 · 10−6 m/s for the matrixmaterial with a range of 2 · 10−6 to 3 · 10−5 m/s. It is also the same orderof magnitude as the mean velocity estimates by Momeyer (2014) (0.7 to1 cm/d at VMC of 0.15).As previously noted, the wet-up takes more time to occur at depth. Thissupports the idea that it takes place in a single permeability vertical flowsystem. The farther down from the moisture source, the more dispersed is39Figure 2.13: C3W2 - Moisture content from TDR probes on tip face 3of the covered pile (C3), 2 m west of the centerline (W2). Dashed linesemphasize multi-year trends.the arrival of the wetting front. Calculated velocities are also consistent withmatrix flow rather than preferential flow and dual permeability. Preferentialflow was observed in uncovered piles under major infiltration events (Neuner,2009). It is not the case in the covered pile especially because of the till layer(better graded material without separate matrix and preferential pathways).If preferential flow was to occur in the Type 1 waste rock layer, it would bedampened out by the till layer, as observed by Marcoline (2008) in a coveredtest pile experiment.40At C2W2, the arrival time of the wetting front at the shallowest probeappears consistent with the temperature of the surrounding material. Thewetting fronts always arrive once the material around the probe is unfrozen(Figure 2.12). Presumably, all the material above up to the surface is alsounfrozen. This pattern is observed at C2E2 from 2008 to 2011. However,from 2013 onwards, data show that the annual wet-up from the top-mostprobe to the deepest one occurs (start to end) before complete thaw of thelayers above. Given the precision of temperature readings and the fact thatthe thermistors are approximately 3 m away from the TDR probes, it couldbe assumed that material around the probes at 5 m depth is slightly warmerthan that near the closest thermistors. It could otherwise indicate that evenif the overlying material is frozen, the material around the TDR probe isthawing and releasing enough water to initiate a wet-up. In that case, waterwould actually come from infiltration events from previous years.Annual Freeze-Thaw Another pattern shown by the TDR data is anannual VMC rapid decline and rise controlled by freeze-thaw. Upon freezing,VMC quickly drops to a very low or minimum value (highlighted in Figures2.11-2.13). The VMC drop is significantly quicker than the one caused bydraining of the pores (desaturation) described earlier. Even though VMCis oscillating a lot once a frozen state is reached, the average remains moreor less constant until it thaws again. VMC values under frozen conditionsare not consistent from one location to the other. It varies from around0.01 up to 0.10. That being said, these values should be taken with somecaution, as the TDR probes have not been designed and tested to measure41Table 2.3: Increase of the annual minimum and annual maximum VMCvalues over the 2008-2011 period. Based on Figures 2.11-2.13.Location Depth Increase of annual minimum VMC Increase of annual maximum VMCC3W2 7.2 m 0.09 (from 0.13 to 0.22) - -10 m 0.09 (from 0.06 to 0.15) - -C2E2 5.1 m 0.03 (from 0.17 to 0.20) <0.015.5 m 0.04 (from 0.18 to 0.22) 0.01 (from 0.25 to 0.26)6.9 m 0.04 (from 0.19 to 0.23) <0.018.7 m 0.09 (from 0.01 to 0.11) 0.11 (from 0.13 to 0.24)C2W2 5.2 m 0.04 (from 0.14 to 0.18) 0.02 (around 0.24)unfrozen VMC at sub-zero temperature. Regardless, the qualitative shapeof the data clearly indicates that the soil and moisture undergo freezingcondition. Figures 2.12 and 2.13 show that freezing periods based on theTDR data are consistent with thermal data.Upon thawing, VMC rises approximately back to the value before freez-ing. It is on average 2.6 times slower than the VMC drop occurring uponfreezing.Multi-year VMC increase During the warming period prior to 2012,some probes show a trend of overall VMC increase. Details are shown inTable 2.3. The increase is obvious at C3W2 because there is no annualwet-up and drainage cycle. Even with an annual cylce at C2E2, the annualmaximum and minimum VMC values also increase from year to year. AtC2W2, a long term VMC increase is only significant at 5.2 m depth.Multi-year VMC decrease As the pile cools down after 2012, most ofthe TDR probes show a trend of overall VMC decrease. Indeed, despitesome occurence of annual freeze-thaw or wet-up and drainage, the annual42maximum VMC (under unfrozen conditions) decreases through time (Fig-ures 2.11-2.13). VMC decreases by roughly 0.01 to 0.03 m3/m3 per year.Other observations− At all locations, VMC never reaches more than 0.25 (Figures 2.11-2.13). This value is consistent with the porosity value of the wasterock matrix (Neuner, 2009).− The annual VMC pattern barely changes at C2E2 top probes (5.1,5.5 and 6.9 m - Figure 2.11). Drainage is much slower than observedelsewhere, as if moisture was pooling on top of a less impermeablestructure. Thermal data show that the underlying material is perma-nently frozen beginning in 2012. It supports the occurrence of such astructure that hinders flow.− For unknown reasons, extremely low VMC conditions are seen priorto 2010 at 7.2 and 9.2 m at C2E2 (Figure 2.11). However, it doesn’taffect the different patterns described earlier.Moisture content in the coverAll the ECH2O probes are well calibrated for temperature but many of them(6 out of 9) show inadequate calibration for moisture content measurements.While the absolute VMC error is up to 0.05 in certain probes, variationpatterns and relative VMC values remain reliable for interpretation. Forinstance, changes in VMC are consistent with the associated temperaturedata. In addition, the difference between minimum and maximum VMC43Figure 2.14: ECH2O probe results at chainage 0+094values do not exceed 0.40, which is a realistic number for porosity of a loosetill.Crest far from heat source As shown on Figure 2.2, three probes arelocated at chainage 0+094. The top probe is in the Type 1 material, some-where above the till layer. The middle one is located within the till, likelyin the upper half of the layer. The bottom probe is also located within thetill layer, but at the base, close to the Type 3 waste rock material. Preciselocation data for these probes is not available.In Type 1 material (Figure 2.14-Type 1), variations of VMC from frozento unfrozen state of about 0.23 are consistent with the Type 1 waste rockmatrix porosity (0.23-0.27) found by (Neuner et al., 2013). When the ma-terial is thawed, it shows considerable daily to weekly oscillations duringthe active season that originates from minor wet-up events associated withprecipitation and infiltration occurring at the surface of the pile.44At the top of till probe (Figure 2.14-Top till), VMC variations fromfrozen to unfrozen state ranges between 0.35 to 0.38. These values are highersince the probe is in the till rather than waste rock. Also, since the probe isdeeper, the small VMC oscillations are dampened out and disappear. Figure2.15 show that upon thawing, VMC peaks at its maximum annual value.Then, VMC remains constant for a certain amount of time that increasesfrom year to year. It remains constant while temperature is still rising,except in 2012. The flat constant VMC plateau lasts for approximately 12,21, 33, 35 and 62 days from year 2012 to 2016 respectively. Afterwards,VMC drops of about 0.09. After a certain time that gets shorter from yearto year (negligible in 2016), VMC drops back to its minimum value as thetill freezes. The pattern shown in Figure 2.15 is an indication of saturatedor near saturated conditions. The increasing duration of such conditionsshow that at this specific location, supplies of water increasingly overcomethe rate at which the material can drain.At the Bottom of till probe (Figure 2.14-Bottom till), temperature rangesbetween -1 to -9◦C. The small VMC variations (0.04 on average from peaksto troughs) basically show unfrozen moisture content moving back and forthon the soil freezing characteristic curve (unfrozen moisture content vs. tem-perature – SFCC), as a result of changing sub-zero temperatures. Hence, itmakes sense that no significant amount of water is flowing there. There isno supply or drainage.Batters As shown in Table 2.2, none of the moisture probes within thebatters experience year round frozen conditions.45Figure 2.15: Unfrozen moisture content in the top of the till layer atSTA 0+094 (Data from Figure 2.14b)Figure 2.16: ECH2O probe results on the East batter, 10 m downslope ofthe crest (E10)46On the east side 10 m downslope (Figure 2.16), data from the probewithin the till layer show two VMC peaks every year upon thawing. Thefirst steep peak occurs in mid July, once temperature reaches the vicinity ofzero degree (Figure 2.17). Thermal conditions are such that free moisturefrom above can flow to the probe. The free water would come from meltedice in layers above and from infiltration events (minor peaks also observed).The first VMC increase is rather sudden. It likely indicates that the amountof ice in the pore space within the till is such that it is not enough toblock water from infiltrating when the temperature is near zero. So, thetill is not saturated at this location. Then, VMC drops as temperaturestill remains around zero (zero curtain effect) and ice melts. VMC likelydecreases because the thawing front and the resulting first wetting front ofthe year are advancing deeper.A month or so later once all the latent heat required is absorbed and allthe ice is melted, temperature starts to rise again. The second VMC peakoccurs at this precise moment. The VMC rise to the second peak is moregradual. It appears to be the maximum value of a slow increase controlledmostly by local ice thaw. Water moving downwards from above could alsobe a component of the increase.Finally, as temperature rises, VMC decreases. This means that the wateris moving downward in comparison to the static condition seen below thecrest (Figure 2.18). Otherwise, VMC varies of up to 0.39 between frozenand unfrozen states, which is consistent with till probes located under thecrest.On the same side of the pile, but at the bottom of the till layer, VMC47Figure 2.17: VMC patterns upon thawing at E10 top of till probedoes not increase much upon thawing. If 2015 data is set aside (peaking over0.6), the annual variation from frozen to unfrozen state is only 0.07. Stillsetting aside 2015, VMC peaks at the same value every year and remains atthis value for an increasing amount of time, the same way as at the top ofthe till at 0+094. Here, the plateau lasts a few hours in 2009 to around 14days in 2014 (Figure 2.19). Assuming that the saturated moisture contentis locally 0.07 to 0.1, that would again be an indicator of moisture build-upabove a frozen barrier located somewhere below that point.On the west side of the pile (Figure 2.20), the limited data that appearsto be reliable at the Top of till probe shows that the material around theprobe is not under saturated conditions. Indeed, a VMC decrease as tem-48Figure 2.18: VMC behaviour within the batter (E10 Top of till) and underthe crest (STA 0+094 Top of till) as temperature risesperature is still rising indicates draining conditions. VMC patterns at thebottom of the till in the west batter are similar to the ones at the top of thetill within the east batter. After the initial VMC increase, it decreases astemperature still rises and until temperature decreases back towards 0◦C.The zero-curtain effect upon freezing increases in the last three years (1month in 2014-2015 to 4.5 months in 2015-2016). The ice content alsoseems to increase towards the end of the data set (0.14 in 2014 to 0.20 in2016). Ice content estimations corresponds to the last VMC measurementimmediately before the soil temperature drops below the freezing point.49Figure 2.19: Bottom of till unfrozen moisture content at E10(10 m downslope on the east batter)Figure 2.20: ECH2O Results in the west batter, 16 and 10 m downslopeof the crest (W16, W10)502.4 Discussion2.4.1 Thermal barrierThe comparison between the covered test pile and the uncovered test pileshows that the cover provides strong insulation and thermal regulation ef-fects, significantly reducing the annual temperature swing at depth. Tem-perature of the waste rock below the till layer is more stable. As noted byPham (2013), temperature fluctuates with much smaller amplitudes thanabove the till layer. Data confirm that it still does so 10 years after theexperiment started.According to Pham (2013), thermal diffusivity of till is slightly greaterthan that of waste rock. Hence, the till layer acts as an insulating layersimply by promoting conduction and shutting down air convection. Themost important implication of the observed insulating effects is a thinneractive layer.2.4.2 Barrier to flow and moisture build upWithin the till below the crestAccording to temperature and ECH2O data around STA 0+094, the bottompart of the till is frozen and remains under no-flow conditions year-round af-ter 2009. There is evidence of a flow barrier located in the cover, within thetill layer. Moreover, moisture behaviour in the upper part of the till showsthat the potential for downward drainage is reduced, which is expected tohappen with the occurrence of a flow barrier located deeper. This inter-51pretation relies on the assumption that thermal conditions do not changesignificantly over years. It is also assumed that there is no internal structuresuch as boulders or local impermeable material that would block downwardinfiltration and thus promote moisture accumulation.The data support the hypothesis of moisture build-up after the onsetof a flux barrier under the crest and within the batters. It is not possibleto state whether all of the water is accumulating, but at least a part ofit is. As a result, volumetric moisture content (VMC) increases towardssaturation. Saturation or near saturation conditions are favourable as itseals the cover to oxygen influx (O’Kane, Januszewski, & Dirom, 2001). Italso potentially constrains the active layer thickness by acting as a latentheat layer (Shiklomanov & Nelson, 2013).Above the previously heated coreFrom 2012 onward, the long term VMC decline in the core above the BCL’sindicates that the infiltration into the core is less than drainage. On onehand, as the whole pile is cooling with the permafrost table rising intothe core, there is no reason to expect an increase in drainage. Indeed, thepile’s outflow falls close to zero in the 2012-2014 period, which confirms thatdrainage is not increasing.On the other hand, the flow regime at the surface is steady throughtime. Precipitation varies annually, but no long term trends over the 10-year period is evident. A decrease in water supply must then come fromflow mechanisms occurring between the topmost TDR probe and the pile’ssurface, namely the cover. A year-round frozen zone would develop within52the cover (at C2W2 and C3W2) and act as a barrier to flow. As a result,water supply to the underlying unfrozen core is hindered. This would cause adecline of VMC over years in material in the vicinity of the moisture probes.Annual wet-up and drainage patterns are not possible anymore.In the case where the frozen zone is not in place year-round, a clearVMC decline over years would not be noticeable. That might be what ishappening at C2E2 location. At this location, it is also possible that someinternal structures prevent drainage. Such a structure could be highlightedby the significantly colder zone previously noted in the thermal results.Below the battersGiven the patterns in ice content variation and VMC upon thawing, it ap-pears that within the active zone in the batters, water is mobile when it isnot under frozen conditions. Peaks in VMC indicate that the extra moistureinvolved upon thawing and infiltration can exit the vicinity of the probe bydraining away.Data from one probe assumed to be closer to the frozen flow barrier isshowing that, just as in the case under the crest, the underlying drainagepotential is reduced. The material remains close to saturation for an in-creasing amount of time year after year, suggesting that moisture build-upalso occurs in slope areas of the pile. This is consistent with the zero-curtaineffect.53Feedback on thermal regimeSignificant moisture build-up at the base of the active layer can form a latentheat layer (Pham, 2013; Stevens, 2016). The amount of latent heat involvedin phase change in such a layer is known to restrain thaw (Carey & Woo,1998; Shiklomanov & Nelson, 2013) or, in some cases, slow down the wholefreeze-thaw process (McFadden & Bennett, 1991; Woo, 2012). It usuallyenhances thermal stability at depth and potentially reduces the active layerthickness (Shiklomanov & Nelson, 2013). The significance of zero-curtaineffects in thermal data does suggest that a latent layer is forming in thecover. However, no effect on the active layer thickness or thawing rates hasyet been observed. Data collected over an extended period would help toassess the extent of the latent layer and its feedback on the active zone.Regardless of any latent heat mechanism, it is likely that the onset ofsaturation conditions seals pore spaces and hence forces heat transfer tooccur by conduction. Hindering any possible air convection across the cover,it is expected to improve the performance of the thermal barrier discussedearlier.2.4.3 Thermal influence of external heating sourceThe influence of the artificial heating source is significant especially in thecore. However, even if the latter still remains unfrozen, a flow barrier appearsto work regardless, as long as a part of the cover stays below 0◦C year-round.This phenomenon, shown mostly by TDR data, could also apply in thecontext of heat generating material. A thermal cover could work provided54that heat generation is assessed and taken into account for the choice ofcover’s materials and thickness.2.4.4 HeterogeneityThe results presented in this chapter address different zones of the coveredpile: the vicinity of 0+094 boreholes, Type 3 core below the crest, Type 3core batters, till layer below the crest, till in the batters, west side, eastside, etc. Data shows an overall picture that is similar within each of thesezones. However in some zones, it also points to some discrepancies betweenlocations that reflect of the pile’s heterogeneity.For instance, within the Type 3 core below the crest, freezing patternsindicate that the scale of heterogeneities is more or less in the range ofdistances between neighbouring measurement points. The cold structureat C2E5 (Figure 2.3) does not appear to be more than a meter thick atthat location. As it is not seen at other faces (4-5m away) nor on the westside of the same face (10m away), the feature is probably only up to a fewmeters wide. This is consistent with the fact that most of moisture responsesmeasured at TDR probe lines are consistent with the closest thermistorstring (3m apart) but show some variation from one location to another(TDR probe lines are 4-5m apart). Also, boulders up to 3.5m wide werefound during deconstruction of an uncovered pile at the site (Barsi, 2017).In light of the above, the size of the test pile and the instrumentationcoverage were large enough with respect to the scale of heterogeneities. Theinterpretations are based on an experiment with a scale sufficiently greaterthan the representative elementary volume of the material in place.552.4.5 Proposed conceptual modelAs of the end of 2016, a flow barrier within the till layer below the crestis formed in the vicinity of chainage 0+094. Evidence of moisture build-upshow that the barrier is effective. A thin and growing flow barrier is takingplace in the area above the BCL, on the west side. The barrier has still yetto develop at other locations such as above the BCL towards the east side(Figure 2.21).In the batters, a flow barrier appeared to be forming but still growingnear the crest and the toe. It occurs in the bottom part of the till layer orthe top part of the core Type 3 material. Although moisture builds up atsome locations, material wet-up and drainage indicate that water is movingwithin the batter active zone. The onset and thickness of the frozen flowbarrier as well as moisture build-up processes are overall more significantwithin the east batter than within the west one.The covered pile features a 3-stage thermal cover that evolves unevenlythroughout the pile, regardless of whether there is an artificial heat sourceor not:1. Thermal regulation and insulation for accelerated cooling2. Onset of a frozen barrier to flow at the base of the active layer. Barrierdepth is variable and the material is not necessarily saturated. It ismost likely unsaturated and somewhat permeable to air.3. Moisture (ice and water) build-up towards saturation. A latent heat56Figure 2.21: Integrated thermal hydrologic conceptual model57layer takes place over the flow barrier. Saturated conditions shouldeventually form a barrier to oxygen transport.In the long run, additional cooling is expected to occur in the next 5-20 years. Where it is already in place, the barrier is expected to thicken.Some thickening would proceed outwards, but most of it is likely to occurinwards, as the permafrost table rises and the core cools down until it finallyfreezes. A frozen barrier to flow is expected to form where it has not formedyet. Moisture build-up below the crest is expected to increase. Moisturebuild-up within the batters might be less substantial because of the slopes.It is uncertain that it could be enough to seal the batters to airflow andoxygen diffusion.Regardless of the size of the pile or the occurrence of heat generationwithin the core material, a proper cover can provide thermal insulation andpromote the onset of a perennially frozen layer below a relatively shallowactive zone. Even if it would take years or decades for the core of a fullscale pile to freeze, the quick onset of a flow barrier would: (1) preventmore water to infiltrate into the potentially acid generating material. Aconsequent moisture build-up that saturates a part of the cover consistentlyover the whole pile (if that is achievable geometry-wise) would in addition(2) prevent airflow and oxygen supply to acid generating reactions alreadyoccurring in the potentially unfrozen core.582.5 Conclusion and recommendationsAcquiring thermal and hydrological data over an extensive period (10+ years)in a remote cold climate environment is challenging for many reasons: coldtemperature effects on instrumentation and calibration; instrumentationlifespan; instrumental drift; installation and maintenance work complicatedby frost; data gaps caused by site access or outages; etc. Fortunately thewaste rock test pile has been well instrumented at sufficient locations toprovide meaningful insights to the long term behaviour of its cover:− Heating due to heat trace at the bottom of the pile had a significanteffect on the thermal and hydrological evolution of the pile, especiallyprior to 2012.− Regardless of the heat effect, the cover plays an important role forthermal regulation of the pile and hence, the onset of a frozen barrierto flow.− Once the barrier is formed, moisture is expected to build up withinthe upper Type 1 waste rock layer.− As of 2016, the data indicate that the active layer is 2.5-3.5 m thickbelow the crest at 0+094. In the batters where a perennially frozenzone developed, the active layer expands deeper, beyond the till layerat some locations. It is thinner near the toe and the crest.In order to get a broader integrated picture of moisture behaviour acrossthe cover and the underlying core, findings of this study should be coupled59in the future with complementary work and data analysis. Similarly as inAnterrieu, Chouteau, and Aubertin (2010) and Pellet, Hilbich, Marmy, andHauck (2016), geophysical characterizations would provide a better insighton the scale of heterogeneity within the pile, moisture content, water geo-chemistry and perhaps ice content distribution. As suggested by Aubertin,Molson, Chapuis, and Cifuentes (2007), a geochemical and oxygen flux re-sponse analysis coupled with temperature and flow would refine the inter-pretations. Finally, numerical simulations of coupled moisture flow and heattransfer could help to better understand the field observations and the un-derlying thermo-hydrological processes. Such numerical investigations arepresented in the following chapter.603Numerical Investigations of Moisture Flowand Heat Transfer in a Covered WasteRock Pile Experiencing Freeze-Thaw3.1 IntroductionModelling of freeze-thaw processes in variably saturated soils is an activefield of research. Many models and underlying theories are being developedand tested in order to address issues and challenges related to cold region en-gineering, artificial freezing, agriculture, climate change and environmentalconcerns (Kurylyk & Watanabe, 2013). In northern Canada, freezing, thaw-ing and permafrost processes are in the foreground of management strategiesfor waste rock produced by mining operations. Risks of acid rock drainagefrom waste rock dumps can be mitigated with multi-layer or insulating cov-ers, depending on site conditions (INAP, 2017) and material availability. Inorder to properly design and build these covers in a cold climate, it is es-sential to have a good understanding of the underlying hydrology and heattransport processes.A large scale experimental waste rock covered pile located in a continuous61permafrost environment in the Northwest Territories was instrumented andmonitored for a period of 10 years, as part of the Diavik Waste Rock Project(Chapter 2, Smith et al. (2013), Neuner et al. (2013), Pham et al. (2013)).The 120 × 80 × 14 m high pile has a potentially acid generating waste rockcore (sulphur content greater than 0.08 wt. %) capped with a 4.5 m thickrock-fill thermal cover. The batters were sloped to a 3H:1V angle. The coverconsists of 1.5 m of till topped by 3 m of run-of-mine waste rock (sulphurcontent less than 0.04 wt. %). Extensive thermal modelling was performedby Pham (2013) regarding the pile and its cover.The objective of this chapter is to provide a better understanding of theevolution of moisture content, fluid flow and temperature in the cover andthe core of the experimental covered pile through one-dimensional numericalmodelling of coupled variably saturated moisture flow and heat transfer. Theaim is to characterize active zone dynamics, ice and moisture accumulation,processes behind the onset of a flow barrier, long term equilibrium withevaporation conditions and feedback of moisture content on the thermalregime. The intent is also to identify the key parameters that have aninfluence on the cover behaviour. Even though many input parameters andboundary conditions are based on field data at Diavik, the goal is not toreproduce field measurements but rather provide a complementary analysisof the cover performance based on numerical investigations.The following sections provide details on numerical implementation andmodel set-up, base case results, a sensitivity analysis, results of long termsimulations accounting for evaporation processes and results of a simulation62without a cover. A discussion includes an integration of the results and anappreciation of the limitations and uncertainties regarding the model.3.2 Modelling3.2.1 Numerical implementationCoupled moisture flow and heat transfer transient simulations in one dimen-sion were performed with the finite element based Soil Vision SVFlux andSVHeat packages (SoilVision Systems Ltd., 2016a, 2016b).Governing equationsFlow The mathematical model for 1D variably saturated flow is given bythe following modified version of Richards’ equation:∂θu∂ψ∂ψ∂t+ ρi∂θiρw∂t= ∂∂z[(Kwγw+ Kvγw)∂uw∂z+Kw +(KTw +KTv)∂T∂z](3.1)On the right hand side, the first term describes liquid water and vapourflow due to pore water pressure (uw) gradient, where Kw and Kv (ms-1) arerespectively the hydraulic conductivity of water and vapour moved by headgradient normalized by the unit weight of water γw (Nm-3). The secondterm describes water flow due to gravity along the vertical z axis. Thethird term describes flow due to a temperature (T ) gradient, where KTw andKTv are hydraulic conductivities (m3s-1K-1) for water and vapour under athermal gradient.The first term on the left hand side accounts for changes in storage where63θu is the volumetric fraction of unfrozen water (liquid and vapour), ψ is thematric suction (kPa) and t is time. Also called specific moisture capacity,∂θu/∂ψ is the slope of the soil-water characteristic curve (SWCC). Thesecond term is related to changes in ice content and can also be expressedas:∂(θ − θu)∂ψ∂ψ∂T∂T∂t= −∂θu∂T∂T∂t(3.2)where ∂θu/∂T is the slope of the soil-freezing characteristic curve (SFCC).This sink term (Eq. 3.2) is responsible for cryosuction. For instance, iftemperature falls, ice content increases (unfrozen water content drops) andconsequently, suction increases. Higher suction may locally amplify porewater pressure gradients and draw water in as a result (cryosuction). Thesuction is then subject to decrease as a result of incoming water. Numericalimplementation of this phenomenon is not simple and convergence may takemore time or fail.To avoid this issue (encountered by the author) the second term on theleft hand side is removed. Ice content is still computed with the SFCC butsuction at sub-zero temperatures remains ice-blind, or not affected by ice-content. Suction is consistent with the SWCC and the total water content(ice and unfrozen) as if it was unfrozen. It is then mostly constant at sub-zero temperatures unless a significant amount of water is entering or exitingthe finite element. The other parameters such as hydraulic conductivity,thermal conductivity and heat capacity remain affected by temperature andcomputed accordingly regardless of the enabled or disabled cryosuction op-tion. Given that hydraulic conductivities of vapour are zero except when64evaporation is applied to a boundary condition, the flow equation used forthe base case scenario becomes:∂θu∂ψ∂ψ∂t= ∂∂z[Kwγw∂uw∂z+Kw +KTw∂T∂z](3.3)Heat transfer The mathematical model used to describe 1D heat transferis:C∂T∂t+(Lf∂θu∂T+ Lv∂θv∂T)∂T∂t− Lv ∂θv∂ψ∂ψ∂t= ∂∂z(λ∂T∂z+ L∗)− (Cwqw + Cvqv)∂T∂z(3.4)On the left hand side, the first term represents changes in energy con-tent (sensible heat), where C is the volumetric heat capacity of the porousmedium (Jm-3K-1). The second term accounts for changes in latent heatof fusion Lf (Jkg-1) of the liquid phase with change in the unfrozen watercontent from change in temperature (SFCC slope), as well as the change inlatent heat of vaporization Lv (Jkg-1) of the vapor phase with change invapor content θv from change in temperature. The third term expresses thechange in latent heat of the vapor phase content with change in soil matricsuction and can be expanded as:− Lv(hrρvsatρv)[∂θw∂ψ+ gωw(n− θw)γwRT]∂ψ∂t(3.5)where hr is the relative humidity, ρvsat the saturation vapor density (kgm-3),65ρv the vapor density (kgm-3), θw the volumetric fraction of liquid water, gthe gravitational acceleration (Nkg-1), n the porosity and R the universalgas constant (kgmol-1K-1).On the right hand side, the first term describes soil sensible heat trans-fer by conduction with λ, the apparent thermal conductivity of the soil(Wm-1K-1), as well as the latent heat transfer by vapor and liquid waterunder thermal, pressure and elevation gradients with L∗ defined asL∗ =(LfKwγw+ LvKvγw)∂uw∂z+(LfKTw + LvKTv)∂T∂z+ LfKw (3.6)The second term on the right hand side of equation 3.4 accounts forthermal convection where Cw and Cv are volumetric heat capacities of theliquid phase and the vapor phase (Jm-3K-1). qw and qv (ms-1) stand forliquid water flux density and vapor flux density.Assuming no evaporation and negligible thermal convection, the heattransfer equation used for the base case scenario becomes:C∂T∂t+(Lf∂θu∂T+ Lv∂θv∂T)∂T∂t− Lv ∂θv∂ψ∂ψ∂t= ∂∂z(λ∂T∂z+ LfKwγw∂uw∂z+ LfKTw∂T∂z+ LfKw)(3.7)Thermodynamic model Coupling of equations 3.3 and 3.7 results in anhighly non-linear system. These equations, as well as the definition of manysoil frozen properties, rely on a common thermodynamic-derived model. TheClausius-Clapeyron Equation (CCE) relates pore water pressure, pore icepressure and temperature in soil at sub-zero temperatures. Many different66versions of the equation exist depending on simplifications and assumptions(Kurylyk & Watanabe, 2013). The one used in the SVHeat code is:ψ = γwLfρwg(T − T0T0)(3.8)where T is the soil temperature and T0 is the freezing point of free water.The CCE is not valid at thermodynamic disequilibrium (Lunardini, 1991;Watanabe, Takeuchi, Osada, & Ibata, 2012). In such disequilibrium, ittends to underestimate unfrozen moisture and hydraulic conductivities infrozen soil (Watanabe & Osada, 2016) especially close to the freezing point.However, disequilibrium should not be a major concern if thawing rates areless than 0.1◦C/h, according to Kurylyk and Watanabe (2013). Hence, theCCE is assumed to be valid.3.2.2 Model GeometryThe model presented in this chapter is a one-dimensional representation ofthe covered waste rock test pile located at the Diavik Waste Rock Projectexperimental site (Smith et al., 2013). The 14-m column model includesa potentially acid generating (PAG) 9.5-m thick waste rock core (referredas Region 3) capped with a cover that consists of a 1.5-m thick loose tilllayer (Region 2) underlying 3 m of non-acid generating (NAG) waste rock(Region 1).3.2.3 Modelling assumptionsThe 1D conceptual model relies on some important simplifying assumptions.67Single permeability Flow is assumed to occur in a single permeabilityand single porosity system. In other words, flow is occurring in the matrixmaterial only, in the two waste rock zones. Preferential flow within macro-pores was observed in uncovered waste rock only under infrequent majorinfiltration events (Neuner, 2009). However, if preferential flow does occurin the upper waste rock layer of the covered pile, it is disrupted by the tilllayer, as shown by field data presented in Chapter 2 and observed by Marco-line (2008) in a similar experiment. The assumption is then mostly relevantto the upper waste rock layer.Cryosuction The effects of moisture migration towards the freezing frontunder pressure gradients caused by freezing are neglected for few reasons.Materials are relatively coarse grained whereas cryosuction is more likely tooccur in low permeability materials (Shiklomanov & Nelson, 2013). Also,the SVHeat cryosuction option could not be used with enough confidenceat the time of modeling. It was a new feature that required further testing(SoilVision System Ltd., personal communication, June 6, 2017). Sincecryosuction is not considered, frost heave or growth of ice lenses are ignoredas well.Vertical flow A vertical flow only column is assumed to be a sufficientrepresentation of the pile anywhere below the crest. Non-vertical moisturemotion can potentially occur along depositional structures within waste rockpiles (Dawood & Aubertin, 2014; Fala, Molson, Aubertin, & Bussie`re, 2005),especially in the batters (Neuner, 2009).68Heat Conduction Conduction is assumed to be the dominant heat trans-fer process. The effect of advective and convective heat transfer is likelynegligible if the magnitude of water fluxes is similar to that of the thaw rate(Lunardini, 1998). Any heat transfer by air movement that would requirea third coupling with air flow equations is ignored. This simplification isconsistent with findings by Pham (2013) and supported by more recent fieldobservations regarding the thermal stabilization effect of the cover (Chap-ter 2).Geochemistry Heat generation in the PAG core due to oxidation of sul-phide minerals is assumed to be negligible. The potentially higher salinityof the pore water is assumed not to influence the freezing point of water.3.2.4 ParametrizationHydraulic propertiesParticle size distributions and statistics on hydraulic properties of wasterock and till are detailed in Neuner (2009), Smith (2009) and Pham (2013).Hydraulic characterization of the waste rock used to build the test piles isprovided by Neuner et al. (2013).SWCC The range of van Genuchten-Mualem parameters found by Neuneret al. (2013) for the waste rock matrix spans approximately one order ofmagnitude variation in suction for a given water content. For modellingpurposes, the parameters were adapted to set up an upper and lower limit forwaste rock SWCC (Figure 3.1a). An average SWCC, within these limits and69based on an in situ wetting-drying experiment was chosen for the base casescenario (Figure 3.1a). SWCC parameters for waste rock found in Hopp,McDonnell, and Condon (2011), Lefebvre, Lamontagne, and Wels (2001),Lefebvre, Lamontagne, Wels, and Robertson (2002) fall within Neuner’srange. Wetting-drying hysteresis effects on SWCC’s were not taken intoaccount in the model.More uncertainty arises regarding the till because no hydraulic character-ization has been conducted. Field moisture measurements (ECH2O probeswithin the till layer) suggest a porosity of around 0.35 to 0.4 (void ratios of0.54-0.6). The void ratios fall into the loose soil category according Viklan-der (1998). Also, the till sampled in Viklander’s experiment contained morefines, but void ratios of most of the looser samples were close to 0.54. Thebase case till SWCC shown on Figure 3.1b is based on the van Genuchtenpredictive model detailed in Chapuis, Masse, Madinier, and Duhaime (2015)with a saturated moisture content of 0.35.Other methods can lead to slightly different curves as illustrated onFigure 3.1b and are tested in the sensitivity analysis presented later. TheMK curve is based on the descriptive Modified Kova´cs (MK) model for non-cohesive soil (Abdelkabir, Bussie`re, Aubertin, & Mbonimpa, 2012; Aubertin,Mbonimpa, Bussie`re, & Chapuis, 2003; Aubertin, Ricard, & Chapuis, 1998;Mbonimpa, Aubertin, & Bussie`re, 2006) and the input MK parameters fromthe equations described in Chapuis et al. (2015). The ‘vG via MK’ curve isa van Genuchten model based on parameters predicted by MK, as detailedby Abdelkabir et al. (2012) and tested by Chapuis et al. (2015).At some locations within the till layer, ECH2O probes provided reliable70(a) Waste rock SWCC(b) Till SWCCFigure 3.1: Soil water characteristic curves (SWCC)71unfrozen volumetric moisture content (UVMC) data under freezing condi-tions. van Genuchten SWCC curves (fit to ECH2O on Figure 3.1) werecomputed to find the best fit of UVMC to temperature for that part of thedata. Curve fitting was done assuming that the CCE is valid and adoptingdifferent porosity values. Interestingly, the curves derived from ECH2O dataare close to the one given by the ROSETTA computer program (Schaap,Leij, & Van Genuchten, 2001) for a material having a similar particle sizedistribution as the actual till.The custom fit curve (on Figure 3.1) is a SWCC shaped with the inten-tion of capturing the predicted SWCC based on particle size (Base casecurve) as well as the high suction behaviour shown by ROSETTA andECH2O data.Hydraulic conductivity The saturated hydraulic conductivity Ksat forwaste rock is set to 9 · 10−6m/s, which is the geometric mean found by Ne-uner et al. (2013) for the matrix material. It is consistent with values used byHopp et al. (2011). For the till, the saturated hydraulic conductivity set to1.4·10−5m/s is based on the Chapuis (2004) method. Permeability functionsfor both till and waste rock are based on Ksat and the SWCC parameters(summarized in Table 3.1) and defined automatically within the code, ac-cording to the van Genuchten-Mualem 1980-1976 estimation method. Thehydraulic conductivity for sub-zero temperatures is discussed later in theThermal Properties section.72Table 3.1: Base Case hydraulic properties (van Genuchten-Mualem pa-rameters)Material Model region θsat α nvG (kPa-1) θres Ksat (m/s)Waste Rock R1, R3 0.25 0.20 2.00 0.040 9 · 10−6Till R2 0.35 0.21 3.74 0.067 1.4 · 10−5Table 3.2: Base case volumetric heat capacity (Pham, 2013)Heat capacity (MJm-3K -1)Material Model region Cfrozen CunfrozenNAG Waste Rock R1 2.3 2.4Tillk R2 2.1 2.5PAG Waste Rock R3 2.1 2.2Thermal PropertiesThermal conductivity Thermal conductivities are based on work byPham (2013) on thermal modeling and investigations conducted at the ex-perimental site. Thermal conductivity of waste rock λWR (Wm-1K-1) isdefined as a function of the degree of saturation Sr using the followingequation:λWR = 0.7 + 1.7( 12Sr1 + 11Sr)(3.9)In the till, it is defined as a frozen thermal conductivity of 3.2 Wm-1K-1and an unfrozen thermal conductivity of 2.9 Wm-1K-1.Heat capacity Volumetric heat capacities are based on previous work byPham (2013). Each material has a set of constant frozen and unfrozen heatcapacities, as shown in Table 3.2.SFCC Many studies report the similarities between SWCC’s and SFCC’s(Bittelli et al., 2003; Black & Tice, 1989; Koopmans & Miller, 1966; Spaans73& Baker, 1996; P. J. Williams, 1964). A common way to define the SFCCis then through empirical relationships that relate the two curves or, simplyvia the CCE, to link temperatures and suction values. The latter approachis used in the base case for both till and waste rock. The phase changetemperature range is set to -0.01 to 0.5◦C for both materials. It is basedon the slope of the SFCC as per recommendations of the SVHeat manual.This setting implies that all phase change and latent heat transfer is forcedto occur within this defined temperature range only.Hydraulic conductivity reduction In order to obtain hydraulic con-ductivity as a function of sub-zero temperatures, permeability functions arecoupled to the CCE the same way SFCC’s are calculated. No further re-duction or impedance factor (Lundin, 1990) is applied, following discussionsand recommendations of Kurylyk and Watanabe (2013) and Watanabe andOsada (2016). The frozen permeability functions are shown in Figure 3.2.Initial conditionsThe initial water content in the column is assigned by layer. Within eachmaterial, the initial water content is assumed to be the same. For the basecase, the initial volumetric moisture content is set to 0.06 for the nodes inthe waste rock, which corresponds to field capacity (Neuner et al., 2013).The degree of saturation is set to 90% in the till, based on field samplescollected during the construction of the cover (Pham, 2013).The initial temperature is constant within the whole column for all sim-ulations. For the base case, the initial temperature of every node is set to74Figure 3.2: Base case frozen hydraulic conductivity of till and waste rock7◦C. This corresponds to the average of the first temperature measurementsfrom thermistors located in the core of the pile.3.2.5 Boundary conditionsThe base case model simulates a total of 3366 days (9.2 years) starting onOctober 14, 2006. The simulated period simply corresponds to the availablefield data for boundary conditions (BC).TemperaturePrescribed temperatures were selected for upper and lower thermal bound-ary conditions, so that the model responds in a manner similar to thatobserved at the field site. At the top of the column, the average daily tem-75Figure 3.3: Base case model boundary conditionsperature TBC-Top is specified according to the following sine function:TBC-Top(t) = −8 + 23sin( 2pi365(t+ 171))(3.10)where t is the number of days after 2006-10-13. The function is a fit tometeorological data from nearby Ekati and Diavik mine weather stations.A sine function is used instead of actual daily temperatures to speed upcomputation by over 250%. A prescribed transient temperature was alsoapplied at the bottom of the column, based on a representative thermistorstring located below the base of the covered test pile. Figure 3.3 shows theboundary conditions used for the base case model.76FlowThe top boundary condition for flow is the a transient specified flux (Neu-mann type boundary condition). It is applied as an atmospheric boundarycondition where evaporation is disabled and daily infiltration is specifiedbased on Penman-Monteith infiltration calculations for 2007 to 2016. Cal-culations are detailed in Zak (2017). To avoid potential convergence issues,infiltration data was modified to make sure that no infiltration occurs whenthe surface temperature is negative. Run off option was enabled to pre-vent potential convergence problems, especially at temperatures close to thefreezing point. This option restricts the pore water pressure at the groundsurface to a maximum of 0 kPa. Any additional water that would cause ahigher pressure is removed of the system.A free drainage boundary condition is applied at the bottom to representthe ability of water to drain out of the pile by the drainage system installedat the base of the test pile (Smith et al., 2013).MeshIn the upper waste rock layer (R1), elements are set to be smaller than 1 cm.In the till (R2), they are set to be smaller than 5 cm. In the core waste rock(R3), the automatic grid generation sets the elements to 7 cm. Near theinterface between two different meshes, the solver automatically computesa gradual element size transition.773.3 Base Case ResultsThis section describes the evolution of the hydrologic and thermal regimesfor the base case simulation.3.3.1 HeatAs shown on Figure 3.4a, the simulated annual temperature swing (ATS)in the cover, from 1 to 4.5 m deep, remains higher than the measured ATS(Chapter 2), but overall the trend is consistent . For instance in the model,the ATS below the cover is significantly less (22%) than the ATS at theground surface. The mean annual temperature profile (Figure 3.4b) in themodel is typical of a conduction dominated system. A small change inslope is seen within the till layer because of a higher thermal conductivity.The average temperatures of the model are up to 3◦C lower than measuredtemperatures (average over the 2008-2016 period). The difference is likelydue to a dimensionality issue. The effect of the batters (approximately90% of the pile’s surface) on the internal temperature is not accounted forin the 1D simulation. Consequently, the absolute temperature values areslightly off, but the effect of material properties and heat transfer processeson relative temperature values is well reflected by the simulated ATS.Given the initial and boundary conditions, it takes a few years for thethermal regime to stabilize throughout the column. The active layer thick-ness (ALT) is stable around 2.8 m in the last 3 years of the simulation, whichfalls in the range of the measured ALT (2.5 to 3.5 m). During the activelayer freeze-back, the ratio of downward freezing to upward freezing (Luo78(a) Annual temperature swing (b) Mean annual temperatureFigure 3.4: Modelled temperature compared to field dataet al., 2014; Woo, 2012) follows an increasing trend. Downward freezingaccounts for 2/3 of the total freezing in year 4 of the simulation, when theALT begins to stabilize and rises up to 80% at the end of the simulation.Similar ratios were observed in the field experiment (Chapter 2).No detailed comparison is made with respect to field data since the goalis not to reproduce field observations but rather get a similar preconditionedthermal regime in order to investigate flow and moisture behaviour in thecover. Also, an artificial heat source had a significant influence on all datasets (Chapter 2), but this artefact is not taken into account in the model.The thermal regime in the column eventually reaches a dynamic equi-librium with respect to boundary conditions, at around 2000 days in theupper part of the column. At all times and all locations in the column, therate of change in temperature never exceeds ±0.025◦C/h. Because it is lessthan 0.1◦C/h, thermodynamic disequilibrium is more likely to be negligible(Kurylyk & Watanabe, 2013) and thus the CCE to be valid.79Figure 3.5: Base case cumulative flux (Negative flux means downwardflow)3.3.2 Flow and moistureGiven the initial contrast in volumetric moisture content (VMC), moisturemoves from the till to the waste rock at the beginning of the simulation. Porewater pressure gradients are such that some moisture even reaches the base ofthe overlying waste rock (WR) layer as shown by the cumulative flux at 3 mdepth (Figure 3.5), but most of the moisture in the till enters the underlyingwaste rock core. As a result, till saturation initially drops whereas saturationin both the upper and core waste rock increases (Figure 3.6).After this initial moisture redistribution, the cover freezes in the firstwinter and saturation remains constant until infiltration occurs during thefollowing active season. Infiltration starts around day 272 in R1 and ap-proximately 90 days later in R2, which indicates that the first wetting frontpropagates downward at a rate of about 3 cm/day. This velocity is con-sistent with field estimates (Chapter 2 & Momeyer, 2014). Once the upper80Figure 3.6: Base case results: saturation through each model regionpart of R3 freezes (around day 125), it remains frozen for the remainderof the simulation. This condition prevents any water supply to the coreand allows the underlying unfrozen material to drain until it freezes entirely(around day 2830) as per the bottom boundary condition.Until approximately day 1200, the active layer is present within thetill layer, allowing surface infiltration to travel into the till. As a result,VMC rises in the till layer during active seasons because the bottom partremains frozen year round. Depending on the magnitude and frequency ofthe infiltration events, the saturation through R1 rises upon the arrival ofthe wetting, but can decline later if the wetting front enters the till layer(see around day 700 on Figure 3.6).After 1200 days, the till remains frozen through its entire thickness forthe rest of the simulation. Then, VMC increases in R1 during each active81(a) Temperature (◦C) (b) Saturation(c) Unfrozen Volumetric Moisture Content (d) Ice Content(e) Total Volumetric Moisture Content (f) Pore-Water Pressure (kPa)Figure 3.7: Base Case results. The black line shows the positition of the0◦C isotherm with depth, as it varies through time - refer to subfigure (a)82season, upon infiltration events. Full saturation (Sr = 100%) is reached atthe base of the upper waste rock around day 2200, which is more or less 6years after the beginning of the simulation. The saturated layer thicknessthen increases at a rate of 0.10 to 0.25 m per year to reach 0.64 m at the endof the 9-year simulation. As shown on Figure 3.7b, the bulk of the saturatedlayer thickness is located within the active layer, which means that moisturein the saturated layer experiences phase change annually. When the activelayer is thawed during the summer, VMC at the base of the upper waste rocklayer reaches 0.25 (full saturation), as shown on Figures 3.7c and 3.7e. Whenthe active layer freezes in the winter, most the water turns into ice (Figure3.7d). The remaining unfrozen moisture content (Figure 3.7c) correspondsto the residual saturation defined by the SFCC, derived from the SWCC.Thus, no air can flow through the saturated layer because the pore spacesare saturated with both ice and residual unfrozen water.The simulation was not long enough to make any conclusions about thefeedback of flow and moisture build-up on the thermal regime. This questionis addressed with long term simulations presented later in this chapter.Figure 3.8 indicates that for 64% of the entire simulation period, fluxeswithin the active layer are negligible (less than 10−14m/s). Higher fluxesreflect the influence of infiltration events, while fluxes greater than 10−7m/sonly occur in the initial stage of the model.Finally, even if the run off option was enabled, no actual run off occurred.83Figure 3.8: Cumulative distribution of maximum daily fluxes3.4 Sensitivity analysisIn this section, the effects of the model element sizes are first examined.Then, the results of a sensitivity analysis on thermal conditions, thermalproperties, flow conditions and hydraulic properties are presented. Theanalysis of the results is based on the comparison of key indicators of thethermo-hydrologic behaviour of the cover. The active layer thickness (ALT)averaged over the last 3 years of the simulation is compared to that of basecase (∆ALT = ALTcase i −ALTbase case). For each region of the model (R1,R2 and R3), the degree of saturation through the region at the end of thesimulation (day 3366) is compared to that of base case (refer to Figure 3.6).In the tables to come, ∆R33661 refers to the difference between the saturationthrough R1 at day 3366 of the case in question and that of the base case. Insome cases, the difference in the initial degree of saturation is also presentedwhen it differs from the base case (noted as ∆R01).The water balance error of all the models presented in this sensitivity84study was always smaller than 1%, and in nearly all the cases, it was lessthan 0.5%.3.4.1 Element sizeFour models with finer meshes in R1 and R2 were tested (Table 3.3). Themesh refinement slightly affects the thermal regime. Indeed, the frozen bar-rier depth varies by a few centimetres depending on the size of the elements.The difference is up to 6 cm between case 5 and the base case. Howevernothing indicates that a finer mesh than applied in the base case providesa more accurate solution for moisture and flow behaviour. In cases 2-5, thedegree of saturation through each region remains within 97.5 and 106% ofthe base case.Table 3.3: Approximate node spacing within each region, ignoring the sizeof elements at layer interfacesNode spacing (m)CASE R1 R2 R3 Total number of nodes1 (Base case) 0.01 0.05 0.07 4732 0.005 0.025 0.07 8033 0.0025 0.01 0.07 14954 0.001 0.005 0.07 34535 0.0005 0.0025 0.07 6765Mesh refinement did not provide typical results where the solution ofmodels with increasingly finer mesh converge towards the solution of themodel with the finest mesh. Discrepancies from one case to another do notappear slowly along the way. Instead, sudden rises or drops occur at specificevents and these punctual errors accumulate over time (Figure 3.9). Theseevents correspond to the passing of thawing and wetting fronts as well as85the initial moisture redistribution. Even though the errors occur under pre-dictable conditions upon thawing, the extent remains somewhat random andunpredictable, regardless of element size, as shown in Figure 3.9. Specificdepths and parameters of Figure 3.9 were chosen to emphasize the incon-sistency between the grid refinement and solution accuracy. This indicatesthat the uncertainty arising from the non-linearity in the system and fromthe numerical solution of the governing equations cannot be overcome bygrid refinement. Interestingly, the coarser grid simulation better performedin terms of water balance (Figure 3.9d).In the light of this comparison, the element size chosen for the base casesimulation is assumed to be adequate for the purpose of this study.3.4.2 Thermal conditionsIn order to assess the importance of thermal conditions on the cover be-haviour, different top and bottom boundary conditions as well as initialtemperatures were tested.Initial conditionsCases with a lower or higher initial temperature than that of the base case(T0) were tested (cases 6, 7 and 8 in Table 3.4). Calculations show thatthe temperature difference between the base case model and theses casesremains always less than 0.5◦C after 1000 days. Within less than 1000days, a dynamic equilibrium between top and bottom boundary conditionssupersedes the initial condition values. As a result, the active later thickness86(a) VMC at -3 m (b) Ice content at -3.5 m(c) Suction (kPa) at -4.5 m (d) Water balance errorFigure 3.9: Element size effects on hydraulic parameters and water balanceat different depths in the column. Refer to Table 3.3 for node spacing. Case1 has the coarser mesh, Case 5 the finer.87and the saturation of the three regions are not significantly affected by theinitial temperature (See cases 6-8 in Table 3.4).Boundary conditionsChanges to boundary conditions are shown in Figure 3.10. The results aresummarized in Table 3.4. Changes in temperature at the top of the pilehave a significant effect on the thermal regime and thus on flow behaviourindirectly. With warmer conditions, the active layer can expand in the tilllayer (Case 11) and into the waste rock below (ALT is greater than 4.5 min case 9). The occurrence and timing of the barrier to flow on which thesaturated layer takes place both depend on the active layer dynamics, whichis mostly controlled by surface temperature. For example, a 5◦C increase inthe temperature at the surface (case 9) leads to 5.08 m thick active layer and,in comparison to base case, a reduction in saturation of the upper waste rocklayer by 0.31 m3 of water per m3 of pore space at the end of the simulation.The results clearly stress the importance of the thermal boundary conditionat the ground surface on the cover performance.The difference between the air and ground surface temperature is usuallyexpressed through N-factors, which are based on annual freezing and thaw-ing indices (Andersland & Ladanyi, 1994; Lunardini, 1981). The base casemodel did not apply N-factors to be consistent with the long-term modelcoupled with evaporation (presented later), because the code do not han-dle N-factors when evaporation is enabled. The N-factors found by Pham(2013) were tested specifically in case 12. The results of case 12 and case 9are similar (Table 3.4). The effect of accounting for the N-factor is observed88(a) Top boundary condition (b) Bottom boundary conditionFigure 3.10: Modified boundary conditions for temperaureTable 3.4: Sensitivity analysis results for thermal conditions, where T0is the initial column temperature, TBC-Top and TBC-Bottom are respectivelythe boundary condition for temperature at the top and the bottom of thecolumn. ALT is highlighted in blue when thinner than the base case. Dif-ferences in saturation are highlighted in light gray when greater than 0.05and in dark gray when greater than 0.1.Case modification ALT ∆ALT ∆R33661 ∆R33662 ∆R33663(m) (m)(m3/m3) (m3/m3) (m3/m3)1 (base case) 2.806 T0 - 4◦C 2.8 0.01 0.01 0 07 T0 - 6◦C 2.81 0.02 0.02 -0.03 0.018 T0 + 3◦C 2.79 -0.01 0 0 09 TBC-Top + 5◦C 5.08 2.28 -0.31 -0.01 0.0910 TBC-Top - (up to 5◦C) 2.08 -0.72 0.12 -0.16 0.0111 TBC-Top + (up to 5◦C) 4.41 1.62 -0.30 0.20 0.0412 N-factors 4.92 2.12 -0.31 0.05 0.0713 TBC-Bottom + 5◦C 3.3 0.5 -0.25 0.33 -0.0114 TBC-Bottom - 5◦C 2.61 -0.19 0.03 -0.02 0.1115 TBC-Bottom = -0.8 ◦C for t>700days 2.75 -0.04 0 0 0.0389to be equivalent to a 5◦C increase of the thermal boundary condition. It issignificant but it does not question the confidence in the base case resultsbecause the latter provided a better and cooler thermal regime solution,closer to field measurements.Applying changes to the bottom thermal condition also affects the coverthermal regime, active layer thickness and moisture behaviour (cases 13-15).However, the effects are less significant than if changes are applied to thetop of the model because the bottom boundary is farther from the cover.The results suggest that heat generation within the core would have thepotential to affect the active layer depth significantly.3.4.3 Thermal PropertiesThermal conductivityAs described in the parametrization section, thermal conductivity and heatcapacity vary with temperature or VMC. In order to test the effect of changesin thermal properties, a multiplying factor is applied to one parameter, onematerial at the time, everything else remaining constant. Factors of 2× and0.5× were used so the values remain realistic. A waste rock thermal conduc-tivity as function of temperature was also tested. Results are summarizedin Table 3.5.The effect of applying a multiplying factor to waste rock thermal con-ductivity (λWR) is more significant than to till thermal conductivity. Forexample, the difference in ALT of cases 18-19 is 8 to 10 times that of cases16-17 (Table 3.5). The magnitude of the effect is likely controlled by the90Table 3.5: Sensitivity analysis results for thermal properties. ALT is high-lighted in blue when thinner than the base case. Differences in saturationare highlighted in light gray when greater than 0.05 and in dark gray whengreater than 0.1.Case modification ALT ∆ALT ∆R33661 ∆R33662 ∆R33663(m) (m)(m3/m3) (m3/m3) (m3/m3)1 (base case) 2.8016 0.5 λTill 2.98 0.19 -0.11 0.14 -0.0117 2 λTill 2.71 -0.09 0.02 -0.03 018 0.5 λWR 1.99 -0.80 0.12 -0.17 -0.0119 2 λWR 4.87 2.07 -0.30 0.13 0.0920 λWR(T) 2.33 -0.46 0.11 -0.15 0.0221 0.5 CTill 2.98 0.19 -0.09 0.12 022 2 CTill 2.59 -0.21 0.04 -0.06 023 0.5 CWR 3.06 0.26 -0.12 0.14 0.0124 2 CWR 2.39 -0.41 0.06 -0.09 -0.01relative layer thickness with respect to the whole system, as well as by thedistance of the material from the daily changing air temperature boundarycondition.In case 20, λWR was defined as a function of temperature instead ofVMC. Based on characterization by Pham (2013), the frozen and unfrozenλWR were set to 3.3 and 1.7 W/mK respectively. Such a parametrizationresults in a higher frozen λWR and a lower unfrozen λWR as shown in Figure3.11. Because it is lower during the active season, the model behaves similarto case 18 (Table 3.5). Thus, the way thermal condition is defined in themodel does influence the results.Heat capacityHigher heat capacities damp the progression of changes in temperature atdepth. It results in a slower moving thawing front during the active season91(a) λWR at 2.995 m depth (b) λWR at 2 m depthFigure 3.11: Waste rock thermal conductivity (λWR) defined as a functionof temperature (Case 20) vs. a function of VMC (Case 1 - Base case)and thus a thinner active layer (cases 22 and 24, Table 3.5). Since the wasterock layer is closer to surface and thicker, waste rock heat capacity changeshave a greater effect on the results than changes in till heat capacity. Theresults are similar to the previous conclusion about thermal conductivity.Overall, the system behaviour is more sensitive to changes in thermal con-ductivity than heat capacity, in accordance with Shiklomanov and Nelson(2013).SFCCRather than compute the SFCC by coupling the CCE and the SWCC, itcan also be defined as more simple linear function. Such SFCC’s are morearbitrary but promote numerical stability during phase change (SoilVisionSystems Ltd., 2016b). Figure 3.12 shows the different SFCC tested. Results92Figure 3.12: Soil-Freezing Characteristic Curves based on SWCC (Basecase) and defined as linear functions.summarized in Table 3.6 show that using the linear SFCC’s leads to a thin-ner active layer. This is likely because latent heat absorption and releaseduring phase change is more gradual, given the slope of the SFCC (Fig-ure 3.12). Also, a lower residual saturation implies that a smaller amountof water will remain unfrozen at sub-zero temperatures. Hence, a greaterproportion of water experiences phase change at every thawing and freez-ing front. Consequently, a larger amount of latent heat transfer makes theactive layer thinner (cases 25-26).93Table 3.6: Sensitivity analysis results for SFCC. ALT is highlighted in bluewhen thinner than the base case. Differences in saturation are highlightedin light gray when greater than 0.05 and in dark gray when greater than0.1.Case modification ALT ∆ALT ∆R33661 ∆R33662 ∆R33663(m) (m)(m3/m3) (m3/m3) (m3/m3)1 (base case) Till θres WR θres 2.8025 Linear WR SFCC 0.01 2.3 -0.50 0 -0.02 0.0126 Linear WR & Till SFFCs 0.001 0.01 2.25 -0.54 0.05 -0.09 0.0127 Linear WR SFCC 2.59 -0.21 -0.08 0.12 0.0228 Linear WR & Till SFFCs 2.48 -0.32 -0.02 -0.04 0.02ConvectionThe base case model has been tested with the convection option enabled(refer to the second term of right hand side of Equation 3.4). In 1D, thisshould be viewed as an advective heat transfer process due to infiltrationthrough the column. No difference was noted in the results. This supportsthe assumption that the magnitude of fluxes is such that sensible heat trans-fer by moisture flow is negligible with respect to temperature conditions andproperties.3.4.4 Flow conditionsSWCC and initial conditionsIn this section, different combinations of SWCC and initial conditions areexamined, based on the uncertainty and the curves presented in Figure 3.1.SFCC’s and permeability functions are also changed indirectly because theyare computed from SWCC’s. Results are presented in Tables 3.7 and 3.8.94Table 3.7: Sensitivity analysis results for SWCC and initial moisture con-ditions. Refer to Figure 3.1 for SWCC’s. ALT is highlighted in blue whenthinner than the base case. Differences in saturation are highlighted in lightgray when greater than 0.05 and in dark gray when greater than 0.1.Case ∆ALT ∆R01 ∆R33661 ∆R02 ∆R33662 ∆R03 ∆R33663Waste Rock Base case29 θres = 0.005 0.04 -0.18 -0.13 -0.08 0 -0.07 -0.13 -0.1230 θres = 0.05 0.04 0.05 0.04 0.01 0 0.03 0.04 0.0431 SWCC : Upper limit n.a. 0.03 0.08 0 -0.07 0.03 -0.0132 SWCC : Lower limit, θSat = 0.16 0.25 n.a. -0.18 -0.02 0 0.03 -0.18 -0.17Till Base case33 θSat = 0.40 0.35 0.01 0 0 0 -0.04 0 034 θSat = 0.25 0.35 0.02 0 0 0.01 0.11 0 035 θSat = 0.20 0.35 0.03 0 0 0.02 0.19 0 036 θres = 0.03 0.067 0.01 0 0 -0.02 -0.10 0 037 θres = 0.01 0.067 0.04 0 0 -0.02 -0.15 0 038 SWCC : vG-M via MK -0.01 0 0 0.09 0 0 039 SWCC : MK 0 0 0.02 -0.01 -0.17 0 0Note that cases 31 and 32 did not reach a steady active layer thickness atthe end of the simulation period.As observed with changes in SFCC, the amount of residual unfrozenwater content within the frozen waste rock (defined by the SFCC derivedfrom the SWCC) has as an effect on the amount of latent heat absorbed orreleased upon phase change. The lower the residual saturation, the thinneris the active layer and vice versa. For example, a residual saturation of 0.005(case 29) results in an active layer 18 cm thinner than that of the base case.Effects of changes in till properties are limited to moisture conditionsin the till without significant effects on the thermal regime, as detailed inTable 3.7 (cases 33-39). The results show that the degree of saturation isgreater in a material with lower porosity values (case 34 and 35). Given95Table 3.8: Sensitivity analysis results for the till SWCC custom fit (Figure3.1 for SWCC’s). Differences in saturation are highlighted in light gray whengreater than 0.05 and in dark gray when greater than 0.1.Till SWCC uWR0 uTill0 ∆u0 θWR0 θTill0 K WR0 K Till0 ∆ALT ∆R01 ∆R33661 ∆R02 ∆R33662 ∆R03 ∆R33663Base case -52.26 -3.05 53.2 0.06 0.32 5 · 10-6 1 · 10040 Custom fit -52.26 -3.05 53.2 0.06 0.32 5 · 10-6 1 · 100 0 0 0 0.02 0 0 041 Custom fit -10.4 -11.2 -0.8 0.13 0.14 5 · 10-3 1 · 10-2 0 0.28 0.02 -0.48 0.11 0.28 0.0142 Custom fit -3.8 -4.8 -1 0.21 0.26 1 · 10-1 4 · 10-1 0.01 0.59 0.02 -0.14 0.12 0.59 0.0143 Custom fit -21.7 -3.39 18.3 0.09 0.31 2 · 10-4 1 · 100 0.01 0.11 0.02 0 0.08 0.11 0the initial moisture redistribution in the column, the greater the residualsaturation, the more water will remain in the till (cases 36 and 37). TheMK curve in case 39 is equivalent to a reduction of the residual saturationvalue (Figure 3.1). The ‘vG-M via MK’ curve in case 38 only affects theinitial till saturation. After the initial moisture redristribution, VMC withinthe till reaches the same residual saturation as in the base case.Results in Table 3.8 illustrate the complexity and interdependence of allthe hydraulic parameters on the evolution of VMC in the till layer. Cumu-lative fluxes at the top and bottom of the till layer (Figure 3.13) indicatethree key processes. First, a lower pore water pressure gradient betweenthe till and waste rock (∆u0 = uTill0 − uWR0 ) facilitates flow from the upperwaste rock (R1) to the till (R2). For example, large initial fluxes from R1to R2 are observed in case 42 (∆u0 = −1) whereas upward initial fluxesare observed in case 40 (∆u0 = 53.2), as shown on Figure 3.13a. Second,more water tends to move from R1 to R2 (Figure 3.13a) if the initial wasterock saturation (θWR0 ) and hydraulic conductivity (KWR0 ) are significantlyhigher, as observed in case 42. Third, a lower initial unsaturated hydraulic96(a) Top of the till (R1-R2 interface) (b) Bottom of the till (R2-R3 interface)Figure 3.13: Effects of SWCC and intial conditions on till layer inflow andoutflow. Negative flux means downward flux. Cases are detailed in Table3.8.conductivity in the till (KTill0 ) makes it more difficult for water to movewithin and out of the till (see case 41 on Figure 3.13b).Saturated Hydraulic conductivityBesides the SWCC parameters, the saturated hydraulic conductivity alsoaffects permeability functions and perhaps the model solution. Table 3.9presents the cases tested and a summary of the results. It indicates thatthe effect of KSat on the active layer thickness is negligible. Like most ofthe other cases, saturation differences occur during the first years of thesimulation, before the onset of a frozen flow barrier in the waste rock layerabove the till. A higher KSat in the till (cases 45-46) reduces the ability oftill material to hold water in the early unfrozen times. With a higher wasterock KSat, drainage of R1 and R3 becomes slightly easier, to the benefit ofVMC in R2 (case 47).97Table 3.9: Sensitivity analysis results for saturated hydraulic conductivity.ALT is highlighted in blue when thinner than the base case. Differences insaturation are highlighted in light gray when greater than 0.05 and in darkgray when greater than 0.1.Case Modified parameter Base case ∆ALT ∆R33661 ∆R33662 ∆R3366344 K TillSat = 1.4 · 10-6m/s 1.4 · 10-5 0.02 0.01 0.04 045 K TillSat = 7.1 · 10-5m/s 1.4 · 10-5 0.02 0 -0.03 046 K TillSat = 1.4 · 10-4m/s 1.4 · 10-5 0 0.01 -0.05 047 K WRSat = 3.0 · 10-5m/s 9 · 10-6 -0.02 -0.01 0.02 -0.0348 K WRSat = 2.0 · 10-6m/s 9 · 10-6 0.02 0.03 -0.06 0.05Overall regardless of KSat, there is no significant difference in the rate ofincrease of the saturated layer thickness. The amount of infiltration is likelythe main driver.Surface infiltrationThe effect of infiltration on flow and moisture behaviour was examined intwo models. The infiltration rate was reduced to 80% in the first one andincreased to 120% in the second one. Such changes in the top boundary con-dition for flow did not affect the ALT. However, it affected moisture build-upas expected. The saturated layer is 38% thicker with more infiltration and38% thinner with less infiltration. High fluxes occur at shallow depth wheninfiltration events are greater, but remain overall within the same order ofmagnitude as the base case. Thus, it is not a concern for advective heattransfer processes.983.5 Long term atmospheric couplingIn all simulations previously presented, evaporation processes were neglectedfor two reasons. First, infiltration estimates were available for the study pe-riod. Second, this simplification significantly reduces computation time. Agiven simulation can take more than 10 times longer to run with atmosphericcoupling. This section provides details on set-up and results of 30-year simu-lations that accounts for evaporation processes. The intent is to characterizethe significance of long term loss of water by evaporation as the upper wasterock wets up above the flow barrier.3.5.1 Numerical implementationThe potential evaporation (PE) is calculated according to Penman (1948)model. The following Wilson-Fredlund-Barbour Limiting Equation (SoilVi-sion Systems Ltd., 2016b) is used to compute actual evaporation (AE):AE = PE(exp(−gωwψ10−fcorrγwRTs)− hr)11− hr (3.11)The empirical factor fcorr is meant to prevent overestimation of AE (Alvena¨s& Jansson, 1997). As recommended by Fredlund, Rahardjo, and Fredlund(2012), a factor of -2 is chosen to account for the coarse grained nature ofwaste rock. The AE calculation invokes the assumption that the soil surfacetemperature is equal to air temperature (SoilVision Systems Ltd., 2016b).99Table 3.10: Climatic parametrizationParameter Field average Reference CalibratedWindspeed 15 km/h Chi et al. (2013) 15 km/hRelative humidity 77% Pham (2013) 90%Net radiation Rn = 13.4 + 83.1 sin( 2pi365 t − 1.3)Pham (2013) 0.425Rn3.5.2 Climatic parametrizationAveraged climatic parameters based on field measurments provided by Chiet al. (2013) and Pham (2013) were found to significantly overestimate AEwhen used in equation 3.11. Hence, relative humidity and net radiationwere modified and calibrated in order to match the simulated infiltration(calculated as the difference between rainfall and actual evaporation) tofield infiltration estimates from water balances on 2 m high lysimeters (Fretz,2013). The modifications summarized in Table 3.10 led to a 0.4% slope errorbetween linear fits to field data and simulation results (Figure 3.14).For the precipitation (PR) input, field data is used for the first 9 yearsof the model. For the following projected 20 years, precipitation data ofa representative year was applied every year. This year was selected sothat: (1) the total annual precipitation is the closest to the mean annualprecipitation over the 9-year dataset and (2) the maximum precipitationevent is the closest to the average annual maximum precipitation.The same sine function used for the base case defines the temperature atthe top of the column. At the bottom, temperature at the end of the basecase simulation (day 3366) is -0.55◦C. For the next 20 years, temperaturelinearly decreases to reach -1◦C at the very end (day 10671).100Figure 3.14: Simulated infiltration fit to data3.5.3 ResultsResults shown in Figure 3.15 indicate that AE increases as moisture avail-ability (Saturation through R1) increases at the top of the column. Thesystem evolves towards an equilibrium state, although it still has not beenreached after 30 years. The AE/PR ratio is getting closer to 1 (Figure 3.15a)as the AE reaches 95 mm/day (Figure 3.15b). The rate of change in AE isdecreasing towards zero, similarly to that of saturation through the upperwaste rock layer (Figure 3.15c).Figure 3.15c shows that saturation through the upper waste rock layerreaches almost 90% after 30 years. The layer is fully saturated from 1 to 3 mdepth according to calculations. It is expected that in the long term, such ahigh moisture content layer would affect the thermal regime and reduce theactive layer thickness, given the amount of latent heat consumed or releasedfor phase change (Pham, Sego, Blowes, Amos, & Smith, 2011; Shiklomanov& Nelson, 2013). However, Figure 3.16 (Base Case air temperature) ratherreveals that active layer thickness is minimal after 9 years and then slowly101(a) Evolution of the ratio of actual evaporation (AE) over precipitation (PR)(b) Evolution of AE (c) Evolution of saturation through R1Figure 3.15: Actual evaporation (AE) results compared to precipitation(PR) and the degree of saturation of the upper waste rock layer (R1)102Figure 3.16: Active layer thickness of the 30-year model accounting forevaporation, with and without air temperature increase (1◦C over the last20 years)increases as moisture starts to build-up within the active layer. The VMCincrease results in a thermal conductivity rise that is great enough to overridethe latent heat effect and thus to thicken the active layer. This phenomenonreported in Shiklomanov and Nelson (2013) is consistent with the significanteffect of thermal conductivity denoted by the sensitivity analysis.3.5.4 Climate warmingA simple 1◦C increase over 20 years of the average temperature has beentested on the same model, everything else remaining constant. The activelayer thickness ends up being 24 cm thicker than the 30-year base casemodel increasing temperature (Figure 3.16). The active layer expands withinthe upper part of the till layer, but the rest of the till remains perenniallyfrozen. This simple test demonstrate that a small increase in temperaturemay become an issue over the long term, especially if the core of the waste103rock pile is not frozen, when heat is generated within the core. There couldbe risks of modest thinning of the perennially frozen layer. An accuratelypredicted rate of surface temperature warming along with a fully coupledheat water flow and air flow would assess more precisely the risk of climatewarming on the cover performance.3.6 No cover simulation3.6.1 Model set-upA cover-free waste rock column was simulated for comparison purposes. Thesame PAG waste rock material was applied to all three regions. Reliabletemperature measurements at 11.2 m depth in the Type 3 uncovered testpile provided data for the bottom boundary conditions (Figure 3.17).Figure 3.17: Bottom temperature condition for the cover-free model fromthermistor data at 11.2 m depth in the T3 uncovered test pile1043.6.2 ResultsThe thermal regime is clearly not accurately simulated, regardless of bound-ary conditions. For example, an active layer thickness of 12 m was observedin the uncovered test pile whereas the model without the cover simulateda 2.6-m thick active layer, which is actually 0.21 m thinner than the basecase. Such a discrepancy is likely a dimensionality issue, as discussed earlierabout the base case results, combined to the fact that airflow is not includedin the model. Indeed, Pham (2013) noted that air convection and wind in-duced cooling occurring below the crest and within the batters significantlycontrol the internal temperatures of the uncovered piles.Otherwise, there is no difference in how moisture behaves within theupper part of the column. Once a frozen barrier to flow takes place at arelatively steady depth, moisture starts building up above it.3.7 Discussion3.7.1 Integration of resultsCover behaviourAccording to the results presented in the previous sections, the cover be-haviour can be characterized in multiple stages:1. During the initial stage, moisture is redistributed within the columntowards a steady state solution. The extent of moisture movementdepends on hydraulic properties, the initial volumetric moisture con-105tent (VMC) in each material and the associated pore water pressuregradients.2. The pile cools down towards an internal dynamic thermal equilibriumcontrolled by temperature variations on the pile surface and at thecontact with the natural ground surface.3. The active layer is thinning until it reaches around 2.8 m and variesby less than 10 cm/year. Unfrozen and frozen properties of waste rockand till are such that the onset of a frozen barrier to flow at the bottomof the active layer is possible.4. A saturated layer forms as moisture builds up on top of the no-flowbarrier. The saturated layer is expected to act as a barrier to oxygendiffusion (O’Kane et al., 2001). Since most of the saturated layer islocated within the active layer, the accumulated moisture is potentiallymobile during the active season, especially if the barrier to flow is nothorizontal. This would not be the case if the till was saturated bycapillary effects before the onset of sub-zero conditions.5. Over the long term, the cover reaches an equilibrium state with atmo-spheric conditions. The net moisture flux across the surface evolvestowards zero as the saturated layer thickness stabilizes. The long termactive layer thickness (ALT) is more affected by a thermal conductiv-ity increase in response to rising VMC rather than by any latent heateffects.106UncertaintyThe model output can vary because of the uncertainty in the input parame-ters such as the SWCC, residual saturation, initial conditions and hydraulicconductivity (range shown by the yellow shading on Figure 3.18). It canalso vary because of the different modelling strategies regarding the thermalconductivity model, the SFCC model and mesh refinement (range shownby the blue shading on Figure 3.18). When to ranges overlap, the shadingappears as green.Figure 3.18 shows that the bulk of the cover behaviour remain the sameoverall. Net moisture flux across the till layer eventually ceases as the frozenmaterial becomes impermeable to flow. The choice of SFCC and thermalconductivity models can delay or move ahead the onset of the flow barrier,which in turns affects saturation in the till. The input parameter uncertain-ties do not affect the thermal regime significantly. Discrepancies in satura-tion are controlled by initial moisture conditions and hydraulic properties.Limitations of the simulation domainOver or underestimation of local moisture build-up might arise from ignoringnon-vertical flow. Estimations from the 1D model are not applicable in zonessuch as slopes, batters or around any non-horizontal depositional structureswithin the cover.It is also possible that the onset of a saturated layer helps to promotethermal insulation by blocking any possible airflow across the cover. How-ever, such a statement could not be verified given that airflow equations107Figure 3.18: Dispersion of the results associated with uncertainty of theinput parameters (yellow) and the modelling approach (in blue)108were not accounted for. It is also expected that the temperatures of the first3 m of waste rock would have been better simulated if wind induced coolingwas accounted for. The latter comment is supported by the results of theno cover simulation.It is important to recall that the heterogeneous nature of waste rock wasnot accounted for in the model. Heterogeneity and preferential flow couldpossibly lead to more variable moisture patterns within the upper wasterock layer under both frozen and unfrozen conditions. However, as far asthe thermal regime remains unaffected, the till layer would likely preservethe integrity of the frozen barrier to flow since it is a more homogeneousmaterial.Key parametersThe rate of moisture build-up and the depth where it takes place are mostlycontrolled by the active layer dynamics, which is strongly dependent on thethermal regime. Sensitivity analysis results clearly indicate that:− The boundary conditions for temperature have the most importantinfluence on the column’s thermal regime. The boundary condition atthe top has a greater control on temperatures in the cover than theone at the bottom.− Within the framework of the aforementioned limitations, thermal con-ductivity is the most influential material property on the cover be-haviour. The effect of thermal conductivity of the waste rock is moreimportant than that of the till.109− Volumetric heat capacity and SFCC also play a significant but lesserrole than thermal conductivity.− Hydraulic properties and conditions do not affect much the thermalregime given the low magnitude of fluxes.− The amount of infiltration and evaporation controls the saturated layerthickness on the long term.− The till can potentially form a latent heat layer that is out of reachof seasonal thaw. The amount of moisture in the till was found to bestrongly controlled by the initial conditions, SWCC’s and permeabil-ity functions. The thermal regime then controls the amount of timeavailable for VMC changes to occur before the till freezes perennially.3.7.2 Implications for cover design and constructionModelling results reported in the present chapter illustrate a few elementsto be considered for design and construction, depending on how the cover isexpected to behave.Barrier to flowIf the goal is to build a thermal cover in which a frozen barrier takes place,cover thickness, thermal properties and N-factors are key elements to con-sider. Since the only possible way to influence thermal properties is bycontrolling moisture content and perhaps soil density (Abu-Hamdeh, 2003;De Vries, 1952; Wierenga, Nielsen, & Hagan, 1969; Yadav & Saxena, 1973),efforts should focus on characterization. The thickness of each layer should110be determined accordingly, so that the active layer does not expand into thematerial underneath the cover.Latent heat layerIf it is aimed to form a perennially frozen high moisture layer within thetill underneath the active zone, the contrast in hydraulic properties betweentill and waste rock should be such that it prevents desaturation of the till.Indeed, modelling showed that even if the till is almost saturated initially,VMC can drop significantly because of a higher initial suction in waste rockand a higher hydraulic conductivity in the till. Compaction could addressthis issue by reducing the hydraulic conductivity of the till so moisturemovement into the underlying waste rock is reduced. The effects of tillcompaction on hydraulic properties are extensively discussed in O’Kane etal. (2001); O’kane, Wilson, and Barbour (1998); Swanson, Barbour, andWilson (1994); Watabe, Leroueil, and Le Bihan (2000). It is also expectedthat compaction would better help to retain water within the sloped tilllayer of the batters. One should be aware that if significant compactionis considered to change cover materials properties, risks of frost heave orice lenses could increase and lead to potential slope failures (Andersland &Ladanyi, 1994; Harris, Smith, Davies, & Rea, 2008; Korshunov, Doroshenko,& Nevzorov, 2016).The cover could also be put in place during the winter. Assuming thatthe till is close to saturation and that the upper part of the waste rockcore is frozen, it could help to keep moisture within the till as long as itremains frozen year round. Otherwise, if the till layer does drain initially,111extra moisture should be provided before it gets frozen year round. Inorder to do so, the active layer needs to penetrate within the till layer for asufficient number of active seasons. To allow this to happen, a proper coverconstruction sequencing should be considered.3.7.3 Comments on freeze-thaw modellingFreezing processes are generally well handled by the code whereas numericalinstabilities were noted regarding thawing. As defined in SoilVision code,phase change is allowed to occur in a range of temperature from Tep (setto -0.5◦C) to Tef (default -0.01◦C) consistently with the SFCC slope. Asa result, the thawing front is in fact a thawing zone (between Tep and Tefisotherms), where many interrelated phenomena are occurring especiallyunder unsaturated conditions.For instance, consider a volume of material within the thawing zone.Initially, a given ice content starts to melt. Numerically, the amount ofunfrozen volumetric moisture content (UVMC) is strictly controlled by theSFCC, which depends on temperature only. But the following mass balancemust always apply :θ = θu + θiwhere θ is the total volumetric moisture content (VMC), θu is the UVMCand θi is the volumetric ice content. At the same time, since warmingoccurs from the top down, the material above the thawing zone is alreadytotally thawed and the water coming from both previously melted ice andinfiltration events is mobile.112The free water can enter the thawing zone below because of unsaturatedconditions (pore space available) and because in theory, unfrozen water canexist at these temperatures. VMC is allowed to increase, but the UVMC stilldepends on temperature, regardless of the incoming moisture. The VMCincrease then results in an ice content increase, as if the incoming waterwould immediately freeze. The ice content increase is especially importantvery close to Tef (Figure 3.19). In some cases, the ice content increasessignificantly above the level it was at before the initiation of thawing (asshown in Figure 3.19a for example). Once Tef is finally reached, all theice is melted and the UVCM (then equal to VMC) is most of the timegreater than the VMC before it froze. Sometimes, VMC in the thawingzone becomes higher than that of the unfrozen material above for a veryshort while. As a result, upward fluxes can be computed.Moreover, the rapid VMC increase that usually occurs upon thawingcan cause a high contrast in suction between the thawing material and thefrozen material at the Tep isotherm. This can result in a computation ofvery large fluxes. Unfortunately, the complex discontinuity created by thedownward progression of the thawing front cannot necessarily be avoided byadapting or refining the mesh of the model, as shown by the element sizesensitivity analysis. Similar problems were also encountered with the freeze-thaw package for Hydrus-1D. Computation time is significantly increasedevery time thawing is involved. Solvers that feature automatic time steppingsystematically reduce the time-step size upon thawing. But depending onthe model set up, convergence might still be extremely difficult or impossibleto achieve when the material is thawing.113(a) Example 1 (Results at 3 m depth)(b) Example 2 (Results at 3 m depth)Figure 3.19: Moisture and ice behaviour vs. temperature upon thawing1143.8 ConclusionThis chapter presents coupled moisture flow and heat transfer long termsimulations of a 1D representation of an experimental waste rock pile locatedin a permafrost environment. The results illustrate the importance of thethermal regime to the active layer thickness, as well as to the onset of asaturated layer and underlying frozen barrier to flow. Thermal propertiesof the cover materials are key elements for cover performance. The surfacetemperature condition is the most important parameter. In the light ofthese findings, cover design must be site specific and account for potentialair temperature increases.Results also show that the initial hydraulic conditions and the hydraulicproperties that drives early stage moisture distributions were such that sat-uration through the frozen till layer remained fairly low (39%). Such con-ditions reduce the potential of this layer to act as high moisture latent heatlayer. Rather, moisture builds up within the active layer and thus expe-riences annual freeze-thaw which makes it potentially mobile during theactive season. Assuming that moisture remains in place, the saturated layerthickness stabilizes over the years and reaches equilibrium with surface evap-oration. Long term simulations indicate that the VMC driven increase inthermal conductivity overrides latent heat effects regarding moisture feed-back on the active layer thickness.Modelling performed well overall, even though some problems were notedregarding thawing dynamics. Two-dimensional modelling could possiblyprovide a more accurate simulation of the thermal regime. It could also lead115to better insights about compaction and the effects of non-horizontal struc-tures such as slopes on lateral moisture flow and accumulation processes.In addition, coupling airflow equations could improve thermal modelling ofthe upper waste rock layer and help to characterize the effect of a saturatedlayer on oxygen diffusion through the pile.1164ConclusionThis thesis provides a detailed characterization of the hydrogeological andthermal evolution of a waste rock test pile (120 m × 80 m × 14 m high)composed of potentially acid generating materials. The test pile, located ina continuous permafrost environment, incorporates a 4.5-m thick rock-fillthermal cover. The findings of this research are based two complemen-tary approaches. First, the covered test pile was instrumented and mon-itored over a 10-year period. Data including meteorological information,soil moisture content, ground temperature, infiltration and outflow rateswas collected, processed and interpreted in order to propose an integratedconceptual model of the pile. Second, numerical investigations were under-taken based on site conditions and material properties. Coupled flow andheat transfer 1D simulations were performed and provided valuable insightsabout the thermo-hydrological behaviour of the cover.4.1 Key findingsThe following list summarizes the key finding of this study.1. The cover effectively acts as a thermal barrier by insulating the core117from seasonal freeze-thaw. Both field monitoring and numerical mod-elling reveal that moisture patterns and active layer dynamics withinthe pile are controlled principally by thermal properties of the mate-rials as well as thermal conditions at the pile surface, within the coreand at the contact with the natural ground surface.2. Turning off the heating cables located at the base of the pile in mid-2011 allowed the permafrost table to rise within the core of the pile.Two years later, flow from the basal drain system stopped as the pilebottom became perennially frozen.3. Field observations supported by numerical simulations indicate thatafter 10 years, the active layer thickness at the top of the pile hasstabilized around 2.5-3.5 m. As of 2016, the active layer expandsdeeper and beyond the till layer in the batters and also below thecrest specifically where the effect of the heating cables was the mostsignificant. However, the active layer thickness at these locations isexpected to lessen because the core of the pile is still cooling down.4. At locations where data show that flow across the cover stops afterthe onset of a variably saturated frozen barrier, there is evidence thatmoisture tends to build-up on top of the barrier. A frozen barrier canform in the cover and block flow even if the core of the pile remain attemperatures above the freezing point. As of 2016, it is not apparentthat the cover frozen barrier has formed uniformly below the crest andwithin the batters, but it is expected to do so in the long term.1185. The initial conditions and hydraulic properties of the till and wasterock are the most important parameters to consider to avoid drainageof moisture from the till to the underlying waste rock in the early stageof the pile. Due to this initial drainage, the till layer in the coveredtest pile may not be functioning as a high moisture content layer.6. Given the semi-arid conditions at the site, long term simulations in-dicate that the hydrology of the cover eventually reaches a dynamicequilibrium where the actual evaporation equals the precipitation, as-suming that lateral flow within the active layer is negligible. Thishappens when the upper waste rock saturated layer is approximately2 m thick.7. The feedback of moisture accumulation in the cover on the subsurfacethermal regime could not be assessed with the 10-year data set. Sim-ulations also show that 10 years is too short to draw conclusions ontis feedback process. However, the 30-year long simulations indicatethat the thermal conductivity rises because of the increasing mois-ture content in the upper waste rock layer. As a result, the activelayer thickens by 1-2 cm a year and stabilizes around 3 m when thedynamic equilibrium in reached in the cover. The effect of thermalconductivity seems to overcome any latent heat layer effects.1194.2 RecommendationsIn the light of the research results and all the observations made throughthe research, the following recommendations are supplied:1. The interpretation of the 1D simulation results is limited by the mod-elling assumptions and also because the pile geometry is dominatedby sloped batters, which cannot be represented by a simple columnmodel. It would be useful to run two-dimensional simulations that ac-count for moisture flow, airflow, heat transfer and climatic conditions.It could provide a better understanding of the cover performance inthe batters and the hydrology of the whole pile. It would help to assesshow mobile is the moisture within the batters and how well the flowbarrier and high moisture layer act as a barrier to oxygen diffusion.2. As mentioned in Chapter 2, it would be interesting to conduct geophys-ical surveys on the test piles. Coupled with punctual measurementsfrom the instruments within the pile, it could help to get a broaderpicture of the cryohydrogeology of the crest and the batters. Variousfield methods of hydrogeophysics combined to inversion models couldprovide valuable insights to the scale of internal structures, porosity,permeability, moisture content, aqueous geochemistry and ground icedistribution.3. Collecting thermal and hydrological data over a decade long periodin a permafrost environment is quite challenging. Cold temperaturesaffect the monitoring instruments, their calibration and lifespan. The120installation and maintenance work is also complicated by frost. All ofthat makes it difficult to acquire continuous and reliable data. Therewould be value in issuing a technical publication compiling the lessonslearned of an exhaustive 10-year waste rock pile monitoring programin a cold region. It would certainly be useful for further researchprograms and professionals that undertake similar work.4.3 Concluding remarksAcid mine drainage is an issue of significant concern when mine wastes con-tain sulphide minerals and limited neutralizing potential. The advancementof research in croyhydrogeology is essential for better predicting and miti-gating the risks of contaminated water discharge to the environment in coldregions. That being said, one should bear in mind that stockpiling wasterock with protective covers is a solution that deals only with the symptomsof a problem to be addressed within a systemic and sustainable perspec-tive. Limiting the production or leaving the ore into the ground should beconsidered by regulatory bodies as viable options for reducing the social,environmental and economic impacts of mine waste production. As men-tioned by a speaker at the 23rd Annual British Columbia-MEND ML/ARDWorkshop, it is important to acknowledge that“not all the ore bodies need to be mined.”121ReferencesAbdelkabir, M., Bussie`re, B., Aubertin, M., & Mbonimpa, M. (2012). Con-version of the modified kova´cs model parameters to the brooks andcorey and van genuchten model parameters for the water retentioncurve of sandy and silty soils. Journal of Irrigation and DrainageEngineering, 139 (5), 388–398.70Abu-Hamdeh, N. H. (2003). Thermal properties of soils as affected bydensity and water content. Biosystems engineering, 86 (1), 97–102.110Akcil, A., & Koldas, S. (2006). Acid mine drainage (amd): causes, treatmentand case studies. Journal of Cleaner Production, 14 (12), 1139–1145.11, 18Al-Houri, Z., Barber, M., Yonge, D., Ullman, J., & Beutel, M. (2009).Impacts of frozen soils on the performance of infiltration treatmentfacilities. Cold Regions Science and Technology, 59 (1), 51–57.10Alvena¨s, G., & Jansson, P.-E. (1997). Model for evaporation, moisture andtemperature of bare soil: calibration and sensitivity analysis. Agricul-tural and Forest Meteorology, 88 (1), 47–56.99Amos, R. T., Blowes, D. W., Bailey, B. L., Sego, D. C., Smith, L., & Ritchie,A. I. M. (2015). Waste-rock hydrogeology and geochemistry. AppliedGeochemistry, 57 , 140–156.11Amos, R. T., Blowes, D. W., Smith, L., & Sego, D. C. (2009). Measurementof wind-induced pressure gradients in a waste rock pile. Vadose ZoneJournal, 8 (4), 953–962.19, 32Andersland, O. B., & Ladanyi, B. (1994). An introduction to frozen groundengineering. Chapman & Hall.2, 5, 9, 10, 25, 31, 88, 111Andersland, O. B., Wiggert, D. C., & Davies, S. H. (1996). Frozen soil sub-surface barriers: Formation and ice erosion. Journal of Contaminant122Hydrology, 23 (1-2), 133–147.9Anterrieu, O., Chouteau, M., & Aubertin, M. (2010). Geophysical charac-terization of the large-scale internal structure of a waste rock pile froma hard rock mine. Bulletin of Engineering Geology and the Environ-ment, 69 (4), 533–548.60Aubertin, M., Mbonimpa, M., Bussie`re, B., & Chapuis, R. (2003). A modelto predict the water retention curve from basic geotechnical properties.Canadian Geotechnical Journal, 40 (6), 1104–1122.70Aubertin, M., Molson, J., Chapuis, R. P., & Cifuentes, E. (2007). Discussionof” variations in moisture content for a soil cover over a 10 year period”.Canadian Geotechnical Journal, 44 (1), 103–106.60Aubertin, M., Ricard, J.-F., & Chapuis, R. P. (1998). A predictive modelfor the water retention curve: application to tailings from hard-rockmines. Canadian Geotechnical Journal, 35 (1), 55–69.70Azapagic, A. (2004). Developing a framework for sustainable developmentindicators for the mining and minerals industry. Journal of CleanerProduction, 12 (6), 639–662.11, 18Azmatch, T. F., Sego, D. C., Arenson, L. U., & Biggar, K. W. (2012). Usingsoil freezing characteristic curve to estimate the hydraulic conductivityfunction of partially frozen soils. Cold Regions Science and Technology,83 , 103–109.3, 4Bailey, B. L., Blowes, D. W., Smith, L., & Sego, D. C. (2015). The diavikwaste rock project: Geochemical and microbiological characterizationof drainage from low-sulfide waste rock: Active zone field experiments.Applied Geochemistry, 62 , 18–34.19Bailey, B. L., Smith, L. J., Blowes, D. W., Ptacek, C. J., Smith, L., & Sego,D. C. (2013). The diavik waste rock project: Persistence of contam-inants from blasting agents in waste rock effluent. Applied Geochem-istry, 36 , 256–270.19Barsi, D. (2017). Spatial variability of particles in waste rock piles (Unpub-lished master’s thesis). University of Alberta.12355Bittelli, M., Flury, M., & Campbell, G. S. (2003). A thermodielectricanalyzer to measure the freezing and moisture characteristic of porousmedia. Water Resources Research, 39 (2).3, 73Black, P. B., & Tice, A. R. (1989). Comparison of soil freezing curve and soilwater curve data for windsor sandy loam. Water Resources Research,25 (10), 2205–2210.3, 73Carey, S. K., & Woo, M.-k. (1998). A case study of active layer thawand its controlling factors. In Proceedings of the Seventh InternationalConference on Permafrost, Yellowknife, Canada (Vol. 10, pp. 23–27).7, 54Chapuis, R. P. (2004). Predicting the saturated hydraulic conductivity ofsand and gravel using effective diameter and void ratio. CanadianGeotechnical Journal, 41 (5), 787–795.72Chapuis, R. P., Masse, I., Madinier, B., & Duhaime, F. (2015). Water-retention curves of coarse soils without organic matter: Improved datafor improved predictions. Geotechnical Testing Journal, 38 (3), 325–337.70Chi, X. (2010). Characterizing low-sulfide instrumented waste-rock piles:image grain-size analysis and wind-induced gas transport (Unpublishedmaster’s thesis). University of Waterloo.20Chi, X., Amos, R. T., Stastna, M., Blowes, D. W., Sego, D. C., & Smith, L.(2013). The diavik waste rock project: implications of wind-inducedgas transport. Applied Geochemistry, 36 , 246–255.19, 100Coulombe, V., Bussie`re, B., Coˆte´, J., & Garneau, P. (2012). Performance ofinsulation covers to control acid mine drainage in cold environments. InCold regions engineering 2012: Sustainable infrastructure developmentin a changing cold environment (pp. 789–799).10Dawood, I., & Aubertin, M. (2014). Effect of dense material layers onunsaturated water flow inside a large waste rock pile: A numericalinvestigation. Mine Water and the Environment, 33 (1), 24–38.68de Grandpre´, I., Fortier, D., & Stephani, E. (2012). Degradation of per-124mafrost beneath a road embankment enhanced by heat advected ingroundwater 1 1 this article is one of a series of papers published inthis cjes special issue on the theme of fundamental and applied re-search on permafrost in canada. Canadian Journal of Earth Sciences,49 (8), 953–962.8De Vries, D. (1952). A nonstationary method for determining thermalconductivity of soil in situ. Soil Science, 73 (2), 83–90.110Elberling, B. (2004). Disposal of mine tailings in continuous permafrostareas: environmental aspects and future control strategies. In Cryosols(pp. 677–698). Springer.10Fala, O., Molson, J., Aubertin, M., & Bussie`re, B. (2005). Numericalmodelling of flow and capillary barrier effects in unsaturated wasterock piles. Mine Water and the Environment, 24 (4), 172–185.68Fredlund, D. G., Rahardjo, H., & Fredlund, M. D. (2012). Unsaturated soilmechanics in engineering practice. John Wiley & Sons.99Fretz, N. M. (2013). Multi-year hydrologic response of experimental waste-rock piles in a cold climate: active-zone development, net infiltration,and fluid flow (Unpublished master’s thesis). University of BritishColumbia.14, 15, 36, 100Gholamzadehabolfazl, A. (2015). A numerical study of a highway embank-ment on degrading permafrost (Unpublished doctoral dissertation).University of Manitoba.10Golder Associates. (1997). Baseline data of climate and surface water hy-drology for the Diavik Diamond Mines EIA (Tech. Rep.). 4260 StillCreek Drive, Burnaby, BC, Canada,V5C 6C6: Golder Associates.13, 20Harris, C., Smith, J. S., Davies, M. C., & Rea, B. (2008). An investigation ofperiglacial slope stability in relation to soil properties based on physicalmodelling in the geotechnical centrifuge. Geomorphology, 93 (3), 437–459.111He, H., Dyck, M. F., Si, B. C., Zhang, T., Lv, J., & Wang, J. (2015). Soilfreezing–thawing characteristics and snowmelt infiltration in cryalfs of125alberta, canada. Geoderma Regional, 5 , 198–208.9Hopp, L., McDonnell, J. J., & Condon, P. (2011). Lateral subsurface flowin a soil cover over waste rock in a humid temperate environment.Vadose Zone Journal, 10 (1), 332–344.70, 72INAP. (2017). Global Acid Rock Drainage Guide (GARD). Retrievedfrom gardguide.com (The International Network for Acid Prevention(INAP))x, 11, 12, 19, 61Iwakun, O., Biggar, K., & Sego, D. (2010). Influence of cyclic freeze-thawon the mobilization of LNAPL and soluble oil in a porous media. ColdRegions Science and Technology, 64 (1), 9–18.9Jacobs, J. A., Lehr, J. H., & Testa, S. M. (2014). Acid mine drainage,rock drainage, and acid sulfate soils: causes, assessment, prediction,prevention, and remediation. John Wiley & Sons.11, 18Jones, J. S. (1981). State-of-the-art report - engineering practice in artificialground freezing. Engineering Geology, 18 (1-4), 313–326.9Koopmans, R. W. R., & Miller, R. (1966). Soil freezing and soil watercharacteristic curves. Soil Science Society of America Journal, 30 (6),680–685.3, 73Korshunov, A. A., Doroshenko, S. P., & Nevzorov, A. L. (2016). The impactof freezing-thawing process on slope stability of earth structure in coldclimate. Procedia Engineering, 143 , 682–688.111Krentz, A. (2014). The hydrologic behaviour of waste rock piles in the cana-dian arctic: snowmelt infiltration and the onset of long term freez-ing in test piles (Unpublished master’s thesis). University of BritishColumbia.24, 33Kurylyk, B. L., MacQuarrie, K. T., & McKenzie, J. M. (2014). Climatechange impacts on groundwater and soil temperatures in cold andtemperate regions: Implications, mathematical theory, and emergingsimulation tools. Earth-Science Reviews, 138 , 313–334.10Kurylyk, B. L., & Watanabe, K. (2013). The mathematical representation of126freezing and thawing processes in variably-saturated, non-deformablesoils. Advances in Water Resources, 60 , 160–177.2, 3, 4, 9, 61, 67, 74, 79Kyhn, C., & Elberling, B. (2001). Frozen cover actions limiting amd frommine waste deposited on land in arctic canada. Cold Regions Scienceand Technology, 32 (2), 133–142.10Lefebvre, R., Lamontagne, A., & Wels, C. (2001). Numerical simulationsof acid drainage in the sugar shack south rock pile, questa mine, newmexico, usa. In Proceedings 2nd Joint IAH-CNC and CGS Ground-water Specialty Conference, 54th Canadian Geotechnical ConferenceSeptember 2001.70Lefebvre, R., Lamontagne, A., Wels, C., & Robertson, A. M. (2002). Ardproduction and water vapor transport at the questa mine. In Proceed-ings of the 9th International Conference on Tailings and Mine Waste,Rotterdam Balkema (pp. 479–488).70Lunardini, V. J. (1981). Heat transfer in cold climates. Van NostrandReinhold Company.88Lunardini, V. J. (1991). Heat transfer with freezing and thawing (Vol. 65).Elsevier.67Lunardini, V. J. (1998). Effect of convective heat transfer on thawing offrozen soil. In Proceedings of the seventh international conference onpermafrost. Canada: Yellowknife (pp. 689–95).69Lundin, L.-C. (1990). Hydraulic properties in an operational modelof frozen soil. Journal of Hydrology, 118 (1), 289 - 310. Re-trieved from http://www.sciencedirect.com/science/article/pii/002216949090264X doi: https://doi.org/10.1016/0022-1694(90)90264-X74Luo, D., Jin, H., Lu¨, L., & Wu, Q. (2014). Spatiotemporal characteristicsof freezing and thawing of the active layer in the source areas of theyellow river (sayr). Chinese Science Bulletin, 59 (24), 3034–3045.5, 6, 9, 27, 78Marcoline, J. R. (2008). Investigations of water and tracer movement incovered and uncovered unsaturated waste rock (Unpublished doctoral127dissertation). University of British Columbia.40, 68Mbonimpa, M., Aubertin, M., & Bussie`re, B. (2006). Predicting the unsat-urated hydraulic conductivity of granular soils from basic geotechnicalproperties using the modified kovacs (mk) model and statistical mod-els. Canadian Geotechnical Journal, 43 (8), 773–787.70McFadden, T. T., & Bennett, F. L. (1991). Construction in cold regions: aguide for planners, engineers, contractors, and managers. John Wiley& Sons.2, 32, 54Momeyer, S. A. (2014). Hydrologic processes in unsaturated waste rock pilesin the canadian subarctic (Unpublished master’s thesis). University ofBritish Columbia.39, 80Natural Resources Canada. (1995). Canada permafrost. (National Atlas ofCanada MCR 4177, 1995; 1 sheet, https://doi.org/10.4095/205314)20Natural Resources Canada. (2008). Annual census of mines, quarries andsandpits. Retrieved from http://www.statcan.gc.ca/pub/16-201-x/16-201-x2012000-eng.pdf10Neuner, M. (2009). Water flow through unsaturated mine waste rock ina region of permafrost (Unpublished master’s thesis). University ofBritish Columbia.14, 22, 39, 40, 43, 68, 69Neuner, M., Smith, L., Blowes, D. W., Sego, D. C., Smith, L. J., Fretz, N., &Gupton, M. (2013). The diavik waste rock project: Water flow throughmine waste rock in a permafrost terrain. Applied Geochemistry, 36 ,222–233.19, 44, 62, 69, 72, 74Nordstrom, D., & Alpers, C. (1999). Geochemistry of acid mine waters.Pages 133 to 156 in: The Environmental Geochemistry of Mineral De-posits, Part A. Processes, Techniques, and Health Issues,(GS Plumlee,MJ Logsdon and LF Filipek, editors). Reviews in Economic GeologyA, 6 .11, 18O’Kane, M., Januszewski, S., & Dirom, G. (2001). Waste rock cover sys-tem field trials at the myra falls operations–a summary of three yearsof performance monitoring. In Proceedings of the 25th annual British128Columbia Mine Reclamation Symposium. Campbell River, BC Techni-cal and Research Committee on Reclamation. BC Ministry of Energy,Mines and Petroleum Resources. September 24th to 27th.52, 106, 111O’kane, M., Wilson, G., & Barbour, S. (1998). Instrumentation and moni-toring of an engineered soil cover system for mine waste rock. CanadianGeotechnical Journal, 35 (5), 828–846.111Osterkamp, T. E., & Romanovsky, V. E. (1997). Freezing of the active layeron the coastal plain of the alaskan arctic. Permafrost and PeriglacialProcesses, 8 (1), 23–44.5, 6, 29Outcalt, S. I., Nelson, F. E., & Hinkel, K. M. (1990). The zero-curtaineffect: Heat and mass transfer across an isothermal region in freezingsoil. Water Resources Research, 26 (7), 1509–1516.7, 31Pellet, C., Hilbich, C., Marmy, A., & Hauck, C. (2016). Soil moisture datafor the validation of permafrost models using direct and indirect mea-surement approaches at three alpine sites. Frontiers in Earth Science,3 , 91.60Pham, H. N. (2013). Heat transfer in waste-rock piles constructed in a con-tinuous permafrost region (Unpublished doctoral dissertation). Uni-versity of Alberta.ix, 10, 12, 16, 51, 54, 62, 69, 73, 74, 88, 91, 100, 105Pham, H. N., Sego, D., Blowes, D., Amos, R., & Smith, L. (2011). Diavikwaste rock project: thermal transport in a covered waste rock testpile..101Pham, H. N., Sego, D. C., Arenson, L. U., Blowes, D. W., Amos, R. T., &Smith, L. (2013). The diavik waste rock project: Measurement of thethermal regime of a waste-rock test pile in a permafrost environment.Applied Geochemistry, 36 , 234–245.14, 19, 22, 62Schaap, M. G., Leij, F. J., & Van Genuchten, M. T. (2001). Rosetta: acomputer program for estimating soil hydraulic parameters with hi-erarchical pedotransfer functions. Journal of hydrology, 251 (3), 163–176.72Shiklomanov, N., & Nelson, F. (2013). Thermokarst and civil infrastructure.129Treatise on geomorphology. Academic, San Diego, 354–373.7, 8, 9, 31, 52, 54, 68, 92, 101, 103Smith, L. J. (2008). Test piles project 2006 construction summary (Tech.Rep.). 38 Sheep River Cove, Okotoks, AB T1S 2G7, Canada: FDAEngineering, Ltd. (submitted to Diavik Diamond Mines, Inc., P.O.Box 2498, 5007 50th Ave., Yellowknife, NT, X1A 2P8, Canada)21Smith, L. J. (2009). Building and characterizing low sulfide instrumentedwaste rock piles: Pile design and construction, particle size and sulfurcharacterization, and initial geochemical response (Unpublished mas-ter’s thesis). University of Waterloo.69Smith, L. J., Moncur, M. C., Neuner, M., Gupton, M., Blowes, D. W.,Smith, L., & Sego, D. C. (2013). The diavik waste rock project:Design, construction, and instrumentation of field-scale experimentalwaste-rock piles. Applied Geochemistry, 36 , 187–199.10, 14, 19, 20, 21, 22, 24, 62, 67, 77SoilVision Systems Ltd. (2016a, December). SVFluxTM1D/2D/3D Sat-urated/Unsaturated Finite Element Groundwater Seepage Modeling,Theory Manual. SoilVision Systems Ltd..63SoilVision Systems Ltd. (2016b, February). SVHeatTM1D/2D/3D FiniteElement Freeze/Thaw Modeling, Theory Manual. SoilVision SystemsLtd..63, 92, 99Spaans, E. J., & Baker, J. M. (1996). The soil freezing characteristic: Itsmeasurement and similarity to the soil moisture characteristic. SoilScience Society of America Journal, 60 (1), 13–19.3, 73Stevens, C. (2016, November). Thermal cover design and performancefor closure of mine waste facilities. Retrieved from http://www.na.srk.com/sites/default/files/file/SRK-ChristopherStevens-Thermal Cover Design and Performance-AMA2016.pdf(SRK North America. Presentation at Alaska Miners AssociationConvention)10, 12, 13, 54Sumgin, M. (1927). Permafrost soils in the USSR. Vladivostok.31Swanson, D. A., Barbour, S. L., & Wilson, G. W. (1994). The modellingof moisture movement in engineered soil covers for mine waste rock130dumps. Retrieved from https://open.library.ubc.ca/cIRcle/collections/59367/items/1.0042259 doi: http://dx.doi.org/10.14288/1.0042259111Tan, X., Chen, W., Wu, G., & Yang, J. (2013). Numerical simulations ofheat transfer with ice–water phase change occurring in porous mediaand application to a cold-region tunnel. Tunnelling and UndergroundSpace Technology, 38 , 170–179.10Urban, F. E., & Clow, G. D. (2017). DOI/GTN-P Climate and Active-Layer Data Acquired in the National Petroleum Reserve - Alaska andthe Arctic National Wildlife Refuge, 1998-2015. US Department of theInterior, US Geological Survey.8Viklander, P. (1998). Permeability and volume changes in till due to cyclicfreeze/thaw. Canadian Geotechnical Journal, 35 (3), 471–477.70Wagner, A. M. (2013). Creation of an artificial frozen barrier using hybridthermosyphons. Cold Regions Science and Technology, 96 , 108–116.9Walvoord, M. A., & Kurylyk, B. L. (2016). Hydrologic impacts of thawingpermafrost - a review. Vadose Zone Journal, 15 (6).10Wang, L., Wang, W., & Yu, F. (2015). Thaw consolidation behavioursof embankments in permafrost regions with periodical temperatureboundaries. Cold Regions Science and Technology, 109 , 70–77.10Watabe, Y., Leroueil, S., & Le Bihan, J.-P. (2000). Influence of com-paction conditions on pore-size distribution and saturated hydraulicconductivity of a glacial till. Canadian Geotechnical Journal, 37 (6),1184–1194.111Watanabe, K., & Osada, Y. (2016). Comparison of hydraulic conductivity infrozen saturated and unfrozen unsaturated soils. Vadose Zone Journal,15 (5).2, 67, 74Watanabe, K., Takeuchi, M., Osada, Y., & Ibata, K. (2012). Micro-chilled-mirror hygrometer for measuring water potential in relatively dry andpartially frozen soils. Soil Science Society of America Journal, 76 (6),1938–1945.13167Weeks, B., & Wilson, G. W. (2005). Variations in moisture content for a soilcover over a 10 year period. Canadian Geotechnical Journal, 42 (6),1615-1630. Retrieved from https://doi.org/10.1139/t05-070 doi:10.1139/t05-07019Weismu¨ller, J., Wollschla¨ger, U., Boike, J., Pan, X., Yu, Q., & Roth, K.(2011). Modeling the thermal dynamics of the active layer at twocontrasting permafrost sites on svalbard and on the tibetan plateau.The Cryosphere, 5 (3), 741–757.8Wierenga, P., Nielsen, D., & Hagan, R. (1969). Thermal properties of a soilbased upon field and laboratory measurements. Soil science society ofAmerica Journal, 33 (3), 354–360.110Williams, D. J., Stolberg, D. J., & Currey, N. A. (2006). Long-term mon-itoring of kidston’s “store/release” cover system over potentially acidforming waste rock piles. Proceedings of the 7th ICARD, St Louis,Mo., USA, 26–30.19Williams, P. J. (1964). Unfrozen water content of frozen soils and soilmoisture suction. Geotechnique, 14 (3), 231–246.3, 74Woo, M.-k. (2012). Permafrost hydrology. Springer Science & BusinessMedia.5, 7, 27, 31, 54, 79Yadav, M., & Saxena, G. (1973). Effect of compaction and moisture contenton specific heat and thermal capacity of soils. Journal of the IndianSociety of Soil Science, 21 (2), 129–132.110Zak, J. (2017). Hydrologic investigations of waste rock test piles in apermafrost environment (Unpublished master’s thesis). University ofBritish Columbia.27, 36, 77Zhang, H., Zhang, J., Zhang, Z., Chen, J., & You, Y. (2016). A consolidationmodel for estimating the settlement of warm permafrost. Computersand Geotechnics, 76 , 43–50.10Zhao, Y., Si, B., He, H., Xu, J., Peth, S., & Horn, R. (2016). Modeling ofcoupled water and heat transfer in freezing and thawing soils, inner132mongolia. Water , 8 (10), 424.9Zona, D., Gioli, B., Commane, R., Lindaas, J., Wofsy, S. C., Miller, C. E.,. . . others (2016). Cold season emissions dominate the arctic tundramethane budget. Proceedings of the National Academy of Sciences,113 (1), 40–45.8133


Citation Scheme:


Citations by CSL (citeproc-js)

Usage Statistics



Customize your widget with the following options, then copy and paste the code below into the HTML of your page to embed this item in your website.
                            <div id="ubcOpenCollectionsWidgetDisplay">
                            <script id="ubcOpenCollectionsWidget"
                            async >
IIIF logo Our image viewer uses the IIIF 2.0 standard. To load this item in other compatible viewers, use this url:


Related Items