Thermal Conductivity of Sand–Silt Mixtures

Heat flow in soils affects the design and management of engineered geosystems, such as thermal solar energy storage facilities (Brosseau et al. 2005), shallow and deep thermal foundations (Laloui et al. 2006), hybrid renewable geothermal systems (Bajpai and Dash 2012; Bidarmaghz and Narsilio 2018), gravity-assisted steam flooding for heavy oil recovery (Wang and Dusseault 2003), and nuclear waste disposal strategies (Madsen 1998; Tang et al. 2008; Gens 2010). These geosystems experience a wide range of stress and temperature conditions in addition to changes in fluid chemistry, pressure, and saturation.

Conduction governs heat flow in soils where small pores ($D_{50}<6\mathrm{mm}$) and moderate temperatures ($T<400°\mathrm{C}–500°\mathrm{C}$) diminish the relative effects of free convection and radiation (Yagi and Kunii 1957; Johansen 1977; Farouki 1981). In fact, the number and so-called quality of contacts determine heat transport in soils; in turn, they reflect the state of stress, packing density, mineralogy, and water content (Salmone and Kovacs 1984; Beziat et al. 1988; Ould-Lahoucine et al. 2002; Weidenfeld et al. 2004; Yun and Santamarina 2008; Wallen et al. 2016).

The coordination number is the number of grain-to-grain contacts per particle in a granular mass. The porosity $n$ and coordination number $cn$ for monosized, spherical grains are bounded by simple cubic ($n=0.476$, $cn=6$) and tetrahedral ($n=0.260$, $cn=12$) packings (Deresiewicz 1958). Binary mixtures give rise to more complex granular packings (McGeary 1961). Thus, grain size distribution and packing configuration affect the thermal conductivity of soils. The compounded effects of texture, density, effective stress, and water content on the thermal conductivity of soil mixtures remain unclear. This study explores the evolution of the thermal conductivity of chemically inert sand–silt mixtures as a function of fines content (FC) and effective stress for both air-dry and water-saturated conditions.

Thermal conductivity $k(\mathrm{W}·{\mathrm{m}}^{-1}·{\mathrm{K}}^{-1})$ is the proportionality constant between the temperature gradient $\nabla T({\mathrm{K}·\mathrm{m}}^{-1})$ and heat flux $q{(\mathrm{W}·\mathrm{m}}^{-2})$ across a medium (Fourier’s law):

In gases, thermal conductivity implies energy transfer through random collisions at the molecular scale (kinetic theory). The net thermal energy flux is proportional to the number of molecules $n_m$ per unit volume ($\mathrm{molecules}·{\mathrm{m}}^{-3}$), the molecular root-mean-square velocity $\overline{\mathrm{C}}$ ($\mathrm{m}·{\mathrm{s}}^{-1}$), the average distance traveled between scattering events or mean free path $\lambda$ (m), and the gradient of the average molecular kinetic energy $\nabla \overline{e}$ (J) (Vincenti and Kruger 1965):where the $\theta$ factor is typically 1:3. The average molecular kinetic energy $\overline{e}$ (J) at a given temperature $T$ (K) is a function of the constant volume specific heat $c_V{(\mathrm{J}·\mathrm{kg}}^{-1}·{\mathrm{K}}^{-1})$ and the molecular mass $M$ (kg) (Vincenti and Kruger 1965):

The thermal conductivity of gases follows from Eqs. (2) and (3):

This expression applies to most liquids as well (Broniarz-Press and Pralat 2009), yet fundamental understanding of thermal conduction in liquid water remains relatively undeveloped (Kell 1972; Bresme and Römer 2013; Pang 2014). Water molecules form hydrogen-bonded clusters in liquid water, where the cluster size depends on temperature and pressure (Lenz and Ojamäe 2009; Pang 2014). Kinetic energy propagates effectively through these clusters; thus, water thermal conductivity is four times higher than the thermal conductivity of other nonmetallic liquids, and it is temperature dependent (Kell 1972).

Thermal conduction in solids depends on the type of molecular bonds. In metals, heat transport arises from the movement of free electrons; their collisions with ionic impurities and phonons in the metallic lattice define the mean free path (Mizutani 2001). For insulators and semiconductors, the thermal kinetic energy is transported through random vibrations in the crystal lattice called phonons, where phonon-phonon and mass-difference scattering prevail over phonon-electron scattering in the bulk material at room temperature (Zou and Balandin 2001). Scattering rates determine the mean free path spectrum (Zou and Balandin 2001; Regner et al. 2013; Freedman et al. 2013); for example, phonons with a mean free path $\lambda >1\mathrm{\mu m}$ contribute to 40% of the conductive heat transport in crystalline silicon at room temperature (Regner et al. 2013).

Heat conduction pathways in soils include grains, pore fluids, and grain-to-grain contacts (Yagi and Kunii 1957; Yun and Santamarina 2008). Typical mineral thermal conductivities are on the order of $k_{\mathrm{mineral}}=1$ to ${10\mathrm{W}·\mathrm{m}}^{-1}·{\mathrm{K}}^{-1}$ (Lide 2010); however, the thermal conductivity of dry soils can be one order of magnitude lower due to porosity and thermal resistance at grain contacts, typically $k_{\mathrm{dry}}{<0.7\mathrm{W}·\mathrm{m}}^{-1}·{\mathrm{K}}^{-1}$. The thermal resistance at grain contacts $R_c{(\mathrm{m}}^2·{\mathrm{K}·\mathrm{W}}^{-1})$ is the ratio of the temperature drop $\mathrm{\Delta }{T}_c$ (K) to the heat flux $q_c{(\mathrm{W}·\mathrm{m}}^{-2})$ across contacting surfaces (Bergman et al. 2011):

As illustrated in Fig. 1, the thermal contact resistance arises from the constriction of heat flow lines through contacting asperities, and from phonon boundary scattering when the diameters of contacting asperities approaches the phonon mean free path length (Zou and Balandin 2001; Prasher and Phelan 2006; Prasher et al. 2007). Clearly, contact resistance makes the bulk thermal conductivity grain-size dependent (Tavman 1996; Weidenfeld et al. 2004). Contact area and the distance between surfaces depend upon the normal forces between grains, the mechanical properties of those grains (i.e., strength and stiffness), and the roughness of the contacting surfaces (Greenwood and Williamson 1966; Mikić 1974; Yovanovich 2005). Contact mechanics causes the soil bulk thermal conductivity to vary with the state of stress according to a Hertzian-like power–law relationship (Cortes et al. 2009; Garrett and Ban 2011; Choo et al. 2012; Roshankhah and Santamarina 2014), and its anisotropy reflects the maximum stress direction and fabric anisotropy (Garrett and Ban 2011; Choo et al. 2012).

On the other hand, pore fluids have considerably lower thermal conductivities, e.g., $k_{\mathrm{water}}{=0.6\mathrm{W}·\mathrm{m}}^{-1}·{\mathrm{K}}^{-1}$ for water and $k_{\mathrm{air}}{=0.024\mathrm{W}·\mathrm{m}}^{-1}·{\mathrm{K}}^{-1}$ for air at room temperature and atmospheric pressure (Lide 2010). Overall, the contact resistance is a function of the thermal properties of solids in contact, the solid–solid contact area, the average distance between neighboring surfaces, and interstitial fluid properties, i.e., simply wetting the interstitial space significantly improves the bulk thermal conductivity (Gori et al. 2001; Wallen et al. 2016).

The effective thermal conductivity of granular and porous materials has undergone extensive investigations since the mid-nineteenth century, and reviews are contained in Woodside and Messmer (1961), Progelhof et al. (1976), Farouki (1981), Abdulagatova et al. (2009), and Zhang and Wang (2017). The Appendix presents a summary of selected thermal conductivity bounds, together with some empirical and analytical models.

The arithmetic and harmonic means of the constituent conductivities define the broadest thermal conductivity bounds, originally proposed in 1912 by Wiener (Dagan 1989). Narrower bounds can be obtained by assuming an isotropic, macroscopically homogeneous medium composed of spheres with inner and outer shells of different conductivities where the ratio of the inner-core and outer-shell radii is preserved for all sphere sizes (Maxwell 1873; Hashin and Shtrikman 1962).

There are multiple empirical relationships for the prediction of the thermal conductivity of dry, saturated, and partially saturated soils (Johansen 1977; Ewen and Thomas 1987; Chen 2008). Often, the geometric mean of the constituent thermal conductivities yields an acceptable first-order approximation of the water-saturated soil bulk thermal conductivity (Johansen 1977; Beck 1976). More recent analytical solutions capture salient effects on thermal conductivity, such as effective stress and Hertzian contact (Weidenfeld et al. 2004), interstitial water and diminished contact resistance (Garrett and Ban 2011), and partially saturated spherical packings with water menisci at contacts (see Appendix) (Haigh 2012). Packing assumptions limit the applicability of these analytical solutions.

The specimens tested in this study are mixtures prepared with uniformly graded quarzitic sand ($D_{50}=300\mathrm{\mu m}$, $e_{\min }=0.53$, $e_{\max }=0.89$) and silica flour (99% quartz, $D_{50}=20\mathrm{\mu m}$, $e_{\min }=0.67$, $e_{\max }=1.51$). Table 1 shows the mass proportions for the 22 sand–silt mixtures. Fines content is the mass fraction of silica flour in these sand–silt mixtures defined as $m_{\mathrm{flour}}/m_{\mathrm{total}}$.

Specimens are tested under zero lateral strain conditions [Fig. 2(a)]. Thermal measurements are taken along the centerline of cylindrical specimens (with $7.2\text{-}\mathrm{cm}$ diameter and $7.6\text{-}\mathrm{cm}$ depth). The vertical effective stress ${\sigma }_z^{\prime }$ in the oedometer cell decreases with depth $z$ owing to side friction from the applied boundary value ${\sigma }_{z0}^{\prime }$ (the Janssen effect) (Sperl 2006; Lovisa and Sivakugan 2015)where the arguments in the exponential function include the interfacial friction coefficient between the specimen and the internal cell wall ${\mu }_0$, Jaky’s lateral Earth pressure coefficient $K_0$, and the depth normalized by the cell radius $z/r$ (Lovisa and Sivakugan 2015). Petroleum jelly and a thin polymer wrap are applied to the walls of the oedometer to reduce the interfacial friction coefficient. Here the stress at middepth is adopted as a representative value, which is about ${\sigma }_z^{\prime }=0.83{\sigma }_{z0}^{\prime }$ for ${\varphi }^{\prime }=27°$, ${\mu }_0\approx 0.18$, and $K_0=1-\sin {\varphi }^{\prime }$ (during loading).

Particle size segregation easily occurs during so-called pouring of the binary granular specimen (Valdes and Evans 2008; Richard et al. 1999). To avoid this effect, the specimen is prepared with a small scoop to carefully place each layer ($0.5\mathrm{cm}$ thick) made of similarly mixed sand–silt composition into the thick-walled oedometer cell. Each layer receives 50 drops of a rod falling from a height of 1 cm (100 g, 1 cm in diameter) to ensure a consistently well-packed specimen. The preparation of water-saturated specimens follows a similar procedure with mixtures scooped under water into the submerged cell to minimize size segregation and air entrapment. The soil layer settles for $5\min$ and then undergoes compaction according to the procedure detailed previously before the placement of the next layer. Finally, the specimen mass and initial height are measured to determine the initial dry density. The accumulation of compaction energy with each successive layer negligibly affects the dry mass density with depth because yield induced by low-energy rodding vanishes within a rod diameter and the strain level readily falls below the volumetric threshold strain in these dense specimens. On the other hand, locked-in horizontal stresses are limited to $K_p{\sigma }_z$ at the time of compaction and are readily overtaken by successive loading.

All tests are vertically stress-controlled. The applied vertical stress ranges from ${\sigma }^{\prime }=41\mathrm{kPa}$ to ${\sigma }^{\prime }=\mathrm{2,642}\mathrm{kPa}$ with a stress increment ratio of 2.0. A LVDT installed close to the center of the loading plate measures the specimen vertical deformation used to compute the evolving dry mass density at each loading step. A $10\text{-}\min$ equilibration period is allowed after each loading step before thermal conductivity measurements.

A thermal needle probe is used to measure thermal conductivity. The thermal needle probe inserted into the centerline of the specimen imposes a constant heat flux along its length for 120 s. The thermocouple at the center of the needle measures the temperature response in time. Fig. 2(b) shows a typical temperature signature recorded for an air-dry specimen. The sediment thermal conductivity $k{(\mathrm{W}·\mathrm{m}}^{-1}·{\mathrm{K}}^{-1})$ is a function of the known input heat per unit length $Q{(\mathrm{W}·\mathrm{m}}^{-1})$ and inversely proportional to the slope of the linear region of the temperature-time data plotted on a semilogarithmic scale (Hooper and Lepper 1950; details in Cortes et al. 2009):

After each heat step, $40\min$ elapses to allow the system to return to thermal equilibrium under ambient conditions before each subsequent measurement. A finite difference simulation of radial heat diffusion confirms that the diameter of the specimen is suffciently large for the $120\text{-}\mathrm{s}$ heating phase so that the boundaries do not violate the basic assumptions that underlie Eq. (7). The calibration of the thermal needle probe uses agar-stabilized water to avoid convection, and it is repeated for three different values of the imposed voltage [$\mathrm{V}=1$, 2, and 4 volts of direct current (VDC)—ASTM D5334 (ASTM 2008)]. Results confirm high measurement reproducibility with an error lower than 2%.

Figs. 3(a and b) present the dry mass density as a function of the applied vertical effective stress and FC for the 22 air-dry and water-saturated specimens. Effective stress and saturating pore fluids have only a minor effect on the dry density. Conversely, the FC strongly impacts the dry mass density, which reaches a maximum at around $\mathrm{FC}=0.3$ [in agreement with the revised soil classification system (RSCS) of Park and Santamarina (2017)]. Note that the coarse sand specimens ($\mathrm{FC}=0$) have a higher dry mass density than that of the fine silt specimens ($\mathrm{FC}=1$) for both air-dry and water-saturated cases; this is due to particle shape, roughness, and electrostatic forces, which hinder the packing of fines (McGeary 1961; Cho et al. 2006).

Fig. 4 shows the thermal conductivity data plotted against the applied vertical effective stress. Data trends highlight the effect of stress, FC, and water saturation. The influence of water saturation is most pronounced; in fact, the presence of water increases the thermal conductivity by more than seven times for specimens at the same effective stress and packing density. Both air-dry and saturated specimens exhibit a quasi-linear trend between thermal conductivity and effective stress in semilogarithmic scale, and mixtures with $\mathrm{FC}=0.4$ have the highest thermal conductivity. A reproducibility study showed that the average error in the determination of $k$ is $\pm 6\%$ ($\pm 2\%–7\%$ for air-dry specimens and $\pm 5\%–8\%$ for water-saturated specimens).

Fig. 5 plots the thermal conductivity against the dry mass density for both air-dry and water-saturated specimens at varying FCs across the range of applied stresses. The air-dry and water-saturated mixtures show a moderate increase in thermal conductivity with dry mass density. The effect of stress on thermal conductivity is clearly seen for each specimen along with the associated change in density for each loading step.

Finally, Figs. 6(a and b) show trends in terms of FC. Fig. 6(a) shows the pronounced effect of FC on dry mass density, with a maximum at around $\mathrm{FC}\approx 0.3$, while the influence of effective stress is minor yet consistent. For clarity, the plotted data correspond to the two extreme stress values of ${\sigma }^{\prime }=41\mathrm{kPa}$ and $\mathrm{2,642}\mathrm{kPa}$. Fig. 6(b) depicts the effect of varying the FC on the thermal conductivity of air-dry and water-saturated mixtures at these two extreme stress values. In all cases, thermal conductivity peaks at $\mathrm{FC}\approx 0.4$. Together, these plots highlight the effect of fines and packing density on air-dry and water-saturated thermal conductivity.

The results of this study demonstrate that water saturation is more important for the thermal conductivity of sand–silt mixtures than effective stress and FC (Figs. 4–6). While the thermal conductivity of water is 24 times that of air, it is much smaller than the thermal conductivity of quartz, $k_{\mathrm{quartz}}=6.8\perp$–$12\parallel \mathrm{W}·{\mathrm{m}}^{-1}·{\mathrm{K}}^{-1}$, which varies depending on whether heat flow is perpendicular or parallel to the crystal $c$-axis. Nonetheless, the presence of water significantly improves heat conduction in the interstitial space at grain contacts because water readily fills interstices and forms a meniscus around contacts in hydrophilic grains (Cho and Santamarina 2001; Yun and Santamarina 2008; Tang et al. 2008). Therefore, even low water saturation in the pendular regime has a marked effect on grain-to-grain heat transfer (Wallen et al. 2016).

The effect of FC on thermal conductivity stems from the corresponding increase in packing density and coordination number, i.e., more conduction paths per unit volume (Yun and Santamarina 2008). The packing density in binary mixtures depends on the individual packing densities of the coarse and fine components, the FC (Stovall et al. 1986), and the relative grain sizes (McGeary 1961). The theoretical maximum packing density occurs at $\mathrm{FC}\approx 0.32$ for simple cubic and $\mathrm{FC}\approx 0.21$ for tetrahedral packing when fines are much smaller than the coarse fraction. On the other hand, the minimum packing densities occur when soils contain only fines or coarse ($\mathrm{FC}=0$ and 1). Similar to packing density, the mean coordination number varies with FC and exhibits minima in both the pure coarse and pure fines cases ($\mathrm{FC}=0$ and 1) and a maximum shortly after the fines fully fill the pores of the coarse grains (maximum coarse-to-fine coordination) (Pinson et al. 1998). A broader range of packing densities and coordination numbers exists when the particles are not spherical but angular and rough (Youd 1973; Cho et al. 2006).

Figs. 6(a and b) show the boundaries of the RSCS where the threshold FCs reflect transitions in the mechanical behavior of soil mixtures (Park and Santamarina 2017). The RSCS framework captures the transition in thermal conductivity with FC. The peaks in density and thermal conductivity for both the air-dry and water-saturated mixtures occur near the transition from SF to F classifications ($\mathrm{FC}\approx 0.3–0.4$). Fig. 7(a) illustrates the heat transmission pathways through a few contacts in monosized soils (note that the overall effect of thermal contact resistance in sands is lower than in silts and results in higher thermal conductivity). In mixtures, Fig. 7(b) shows that heat transport takes place through both fine and coarse grains. This effect is most pronounced at $\mathrm{FC}=0.4$ when both fractions are load bearing.

Effective stress improves thermal conductivity through fabric compaction, an increase in coordination number, Hertzian deformation, and the reduction of thermal contact resistance (Garrett and Ban 2011). The dry mass density and the thermal conductivity follow linear relationships with the logarithm of the effective stress (Figs. 3 and 4) (Roshankhah and Santamarina 2014):where ${\rho }_{100}(\mathrm{kg}·{\mathrm{m}}^{-3})$ and $k_{100}(\mathrm{W}·{\mathrm{m}}^{-1}·{\mathrm{K}}^{-1})$ are respectively the dry mass density and thermal conductivity at the vertical effective stress ${\sigma }^{\prime }=100\mathrm{kPa}$. The dimensionless factors $\chi$ and $\beta$ capture the increase in dry mass density and thermal conductivity for a 10-fold increase in the vertical effective stress. Table 1 summarizes the fitted values for all air-dry and water-saturated mixtures.

The thermal conductivity and dry mass density of all mixtures exhibit residual effects upon unloading. The underlying processes include irreversible particle rearrangement, increased coordination number, plastic deformation at contacts, locked-in lateral stresses, and some grain crushing as the effective stress approached yield conditions during loading (Yamamuro et al. 1996; Lade et al. 1996; Sperl 2006).

Data generated during the study are available from the corresponding author upon request.

Support for this research was provided by the Goizueta Foundation at the Georgia Institute of Technology and the KAUST endowment at King Abdullah University of Science and Technology. The authors’ gratitude extends to Gabrielle E. Abelskamp, who edited the manuscript.