Estimation of thermal conductivity of cemented sands using thermal network models

Wenbin Fei,Guillermo A.Narsilio

Department of Infrastructure Engineering,The University of Melbourne,Parkville,Australia

Keywords:Network Cementation Computed tomography Ground improvement Sands

ABSTRACT Effective thermal conductivity of soils can be enhanced to achieve higher ef ficiencies in the operation of shallow geothermal systems.Soil cementation is a ground improvement technique that can increase the interparticle contact area,leading to a high effective thermal conductivity.However,cementation may occur at different locations in the soil matrix,i.e.interparticle contacts,evenly or unevenly around particles,in the pore space or a combination of these.The topology of cementation at the particle scale and its in fluence on soil response have not been studied in detail to date.Additionally,soils are made of particles with different shapes,but the impact of particle shape on the cementation and the resulting change of effective thermal conductivity require further research.In this work,three kinds of sands with different particle shapes were selected and cementation was formed either evenly around the particles,or along the direction parallel or perpendicular to that of heat transfer.The effective thermal conductivity of each sample was computed using a thermal conductance network model.Results show that dry sand with more irregular particle shape and cemented along the heat transfer direction will lead to a more ef ficient thermal enhancement of the soil,i.e.a comparatively higher soil effective thermal conductivity.©2022 Institute of Rock and Soil Mechanics,Chinese Academy of Sciences.Production and hosting by Elsevier B.V.This is an open access article under the CC BY-NC-ND license(http://creativecommons.org/licenses/by-nc-nd/4.0/).

1.Introduction

Effective thermal conductivity(λeff)is one of the key soil parameters required for shallow geothermal engineering.Moisture content is the main factor in fluencing the soil thermal conductivity since water thermal conductivity(0.596 W/(m K))is more than 20 times that of air(0.025 W/(m K))(Young et al.,1996).The thermal conductivity of the solid particles is one order of magnitude higher than that of water.Therefore,soil porosity or bulk density,indicating the fraction of solid particles,also has an in fluence onλeff.Unless the soils are completely saturated(high moisture content),or work in arid environments(low moisture content),their moisture content typically varies and is hard to be controlled;and porosity or bulk density can only be typically enhanced in engineered fills or(typically)in up to very shallow depths(＜1 m).Therefore,innovative approaches are required to enhance the thermal conductivity of soils.Deep mixing and grouting(soil cementation)may constitute one such approach.To understand the in fluence of cementation in heat transfer,and to reduce the number of in fluencing factors(e.g.moisture and compaction),here we study dry geomaterials as a first step since soils in arid environments could be considered mostly dry from their heat transfer mechanisms point of view.

Heat transfer in dry soils at a temperature lower than 700°C(El Shamy et al.,2013;Asakuma et al.,2016)is mainly due to conduction through three paths:particles themselves,interparticle contacts and near-contacts(Yun and Santamarina,2008;Fei et al.,2019a).Grouting,a traditional ground improvement method,can cement pre-existing particles together to achieve a superior bearing capacity for foundations.It also increases the interparticle contact area and thus offers the extra bene fit of superior heat transfer for the ground since interparticle contacts are a primary heat transfer path(Yun and Santamarina,2008;Narsilio et al.,2010a,b).However,the enhancement of soil thermal conductivity for a more ef ficient geothermal system has not been paid enough attention to.Cementation can be achieved by injecting cement grout into the ground(deep mixing and grouting)or by microbially induced calcite precipitation(MICP)technique as an alternative(Wang et al.,2020).The two different approaches may lead to cementation at different locations within the soil matrix or to cementation anisotropy(preferential cementation growth directions),resulting in differentλeffvalues.However,related studies have not been conducted.Additionally,since particles in natural sands have distinct shapes,the investigation of how particle shape affects the enhancement ofλeffdue to cementation is also necessary.

One of the essential works of studying heat transfer in the soils is to calculate or measure theλeffaccurately and ef ficiently.Although laboratory results are reliable,testing requires relatively large non-destructive samples that are sometimes dif ficult or costly to obtain,as for sandy soils in particular.It is more challenging to build an empirical equation for predicting the thermal conductivity of these sandy soils since more samples are required.Another limitation in the empirical derivation is that the microstructure information controlling the macroscopic(or engineering)heat transfer is limited to parameters such as mean particle size,void ratio,and mean coordination number.However,the inaccessibility of local microstructural information such as interparticle contact area and its distribution constrains the fundamental understanding of the mechanisms of heat transfer in natural soils(cemented or uncemented).Even with access to local microstructural information,this information andλefffrom numerous samples are required to build a reliable arti ficial intelligence(AI)predictive model(Zhang et al.,2020;Fei et al.,2021a).

The discrete element method(DEM)has been used to simulate heat transfer in spherical granular media with the convenience of extracting the particle size and coordination number of each particle(Chaudhuri et al.,2006;Moscardini et al.,2018).Nevertheless,it may not be applicable for natural soils in which soil particles show complex shapes even though spheres can be clumped together to imitate an irregular particle.To address this limitation,an X-ray computed tomography(XCT)based level-set DEM(de Macedo et al.,2018)has been developed,but the current version has no functionality of simulating heat transfer.Alternatively,XCT based finite element simulation discretises irregular particles into meshes,in which energy balance equations are solved(Narsilio et al.,2010b;Fei et al.,2019b).Even though the finite element method(FEM)simulates local heat transfer inside the particle and at interparticle contacts,without further XCT image postprocessing during the simulation,FEM cannot deal with particle volume effect and may lead to an overestimation ofλeff(Persson et al.,2005;Wiebicke et al.,2017).

Recently,network models have been introduced to simulate transport properties directly using XCT images but widely used for fluid flow rather than heat transfer(Gostick et al.,2016;van der Linden et al.,2016).A network model presents a sample as a web of nodes and edges.In a pore network,a node presents a void while an edge indicates a pore throat that connects two neighbouring pores.Our research group has pioneered the development of an XCT image based thermal conductance network model(TCNM)that treats particles as nodes linked by edges presenting interparticle contacts or near-contacts(Fei et al.,2019a;van der Linden et al.,2021).TCNM first introduces a greyscale intensity penalty factor for mitigating the overestimation of contact area andλeffdue to the partial volume effect.It also shows high computational ef ficiency by detecting the adjacency of solid and void phases after watershed segmentation,which is hinted at by the dual graph of Voronoi diagram and Delaunay triangulation(Aurenhammer et al.,2013).The TCNM approach and associated tools such as three-dimensional(3D)particle shape descriptor and complex network theory(i.e.graph theory)can achieve richer valuable microstructural parameters(Fei and Narsilio,2020).Hence,this combined approach resolves the hypothesis that thermal front propagation does not merely depend on global porosity(Yun and Evans,2010).The available new network features are important microstructural parameters that can boost more reliable heat transfer models in either a theoretical or AI approach.

This paper aims to investigate the extra bene fit of cementation to enhance not only the geomechanical response of the ground,but also the thermal conductivity of soils.Three natural sands were selected for scanning and cementationwas simulated by dilating the solid phase in XCT images along different directions.Theλeffof the cemented sands was computed using TCNM.With these results,the impact of cementation onλefffor the three sands were compared followed by the impact of cementation anisotropy onλeff.

2.Materials

2.1.Uncemented sands

Glass beads,Ottawa sand and angular sand were selected in this work since the complexity of their particle shape increases from spherical,near-spherical to irregular,as shown in Fig.1.Glass beads and Ottawa sand are two common research materials used in geotechnical experiments.Glass beads are made of silica,while both Ottawa sand and angular sand mainly consist of quartz.Their equivalent particle sizes(Table 1)derived from particle volume are similar:0.6 mm for glass beads,0.76 mm for Ottawa sand,and 0.68 mm for angular sand.XCT images can be freely accessed by readers in Fei et al.(2021b).

Table 1Particle size(mm)of the studied granular materials.

2.2.Cemented sands

Fig.1.The computed tomography(CT)images of three sands:(a)Glass beads,(b)Ottawa sand,and(c)Angular sand.

Grouting is an arti ficial surface nucleating procedure to introduce new solids in the pore space of soils.Natural nucleation is also widely observed in deep geothermal projects(power generation)when steam and fluid are circulated in an open aquifer where mineral precipitation induced by variations of temperature and pressure occurs.To imitate the surface nucleation that generates new solids on the surface of pre-existing particles,a cubic subsample(350 pixels in length,and pixel size is 13μm)from the XCT images of each sand was first cropped,as shown in Fig.2a.It is noticeable that XCT images are two-dimensional(2D)images in sequence with an interval of a pixel size and stacking them can reconstruct a digital 3D sample.Next,Otsu threshold segmentation(Otsu,1979)was utilised to identify the solid and void phases,resulting in a binary image with solid in black and void in white,as shown in Fig.2b.After that,3D dilation was applied to the preexisting solid phase in the binary image to generate the extra solids around the pre-existing particles,as shown in Fig.2c.To better visualise the nucleation,the new solids were rendered as white around the original greyscale particles in Fig.2d.This‘virtual’dilation was varied from 1 pixel to 10 pixels in this work,and the obtained 10 cemented samples of varying cementation degree for each sand were labelled as D1-D10,respectively.

Noticeably,solids could be over-generated during threshold segmentation.The ampli fication of the local images in the yellow rectangles from Fig.2a and b is presented in Fig.3.While solid voxels in Fig.3a are in light grey and real void voxels are in black,some void voxels at the interparticle contacts are partially filled with grey during XCT scanning.These partially filled void voxels were incorrectly grouped as solid in the binary image after threshold segmentation,as shown in Fig.3b,which is called partial volume effect(Wiebicke et al.,2017).In the present work,this overestimated solid phase at the interparticle contact was considered as local cementation(sample D0)to imitate the cementation bridge between soil particles due to mineral precipitation(Wang et al.,2020).Hence,11 cemented samples for each kind of sand were generated to simulate surface cementation in this paper(D0-D10).

MICP is an alternative way to achieve cemented sands(Wang et al.,2020).The underground fluid flow could control the chemical alteration (cementation/mineralisation or dissolution)in different approaches(Lichtner,1988).For example,the propagation of reaction fronts can lead to a banded or layered pattern such as marble.The variations of temperature and pressure along the streamlines could result in gradient reactions.These cases can result in anisotropic cementation along different directions as shown in the yellow windows in Fig.4.The two particles in window A are connected due to horizontal mineral precipitation(Fig.4a),while they are still separated because of vertical mineral precipitation(Fig.4b).In contrast,the two particles in window B keep the small gap after horizontal mineral precipitation(Fig.4a),while they become in contact resulting from vertical mineral precipitation(Fig.4b).The directional cementation can alert the preferential heat transfer path and intensify the thermal anisotropy of the sands.The procedures of producing the samples cemented in all directions have been illustrated in Fig.2.A similar approach but in horizontal and vertical directions was applied to the original sands to generate directionally cemented samples.One-pixel directional cementation was iterated 10 times for each kind of sand.Hence,10 horizontally cemented samples and 10 vertically cemented samples were generated for each kind of sand.

Fig.2.Approach to achieve surface cemented Ottawa sand.The original greyscale CT image(a)was threshold segmented into a binary image(b)followed by a 2-pixel dilation of the solid phase(c),which results in the extra white pixels in(d)presenting grout.

Fig.3.Overestimation of the contact area due to partial volume effect:(a)Original greyscale image;and(b)Binary image after threshold segmentation.

Fig.4.Mineral precipitation in(a)horizontal and(b)vertical directions.A 2-pixel dilation was applied.

3.Thermal conductance network model

3.1.Network construction

The overestimated contact area shown in Fig.3b leads to the overestimation ofλeffsince interparticle contact is a primary heat transfer path,especially in dry sands(Yun and Santamarina,2008).The overestimation ofλeffcannot be mitigated using traditionally simulation methods such as the FEM.To address this barrier,a TCNM was developed by our research team(Fei et al.,2019a;van der Linden et al.,2021)and implemented in the simulation of heat transfer in this paper.TCNM is an XCT and network based simulation tool,thus it provides the bene fit of post-processing images during the computation.

Fig.5.Procedures of simulating heat transfer in the thermal conductance network model(TCNM).

An XCT image stack(Fig.5a)was used to reconstruct 3D digital sand shown in Fig.5b after Otsu threshold segmentation.Anetwork is a web of nodes and edges,and in a thermal network,nodes were used to present particles while edges were used for interparticle contact or near-contacts to cover the main heat transfer paths in dry granular materials(a contact network typically only considers real contacts as edges,not near-contacts).Consequently,watershed segmentation was employed to split the integrated solid phase in Fig.5b into individual particles with unique identi fiers(IDs)and rendered with randomcolours in Fig.5c.For each particle,a black node was assigned to its coordinate centre,as shown in Fig.5d.To detect interparticle contact,boundary voxels of each particle were recognised.Next,if the boundary voxels of a particleiwere shared with particlej,they belonged to interparticle contacts and were presented as red edges,as shown in Fig.5e.To identify near-contacts ef ficiently,watershed segmentation was also applied to the void phase as shown in greyscale in Fig.5d.After that,if particle boundary voxels bordered on a pore,particle ID and pore ID were grouped as a particle pore connection illustrated as a purple arrow in Fig.5d.Then,potential near-contacts shown as blue edge could be ef ficiently identi fied.The near-contacts are essentially a group of gaps between boundary voxels of two neighbouring particles shown as blue cylinders in Fig.5e.If the length of the shortest blue cylinders was smaller than a cut-off rangeε,a near-contact was built in the thermal network as blue edges.In addition,all blue cylinders with lengths smaller thanε will be later used to compute the thermal conductance of the nearcontact.

By now,the thermal network has only presented the topology inside sand.To simulate heat transfer using the thermal network,particles and contacts were treated as conductors connecting the top and bottom boundaries of the whole assembly,as shown in Fig.5f.Therefore,Fourier’s law can be used to calculate the required heat flux transferring through the sample:

whereQijis the heat flux transferring from nodeito nodej,which is induced by a temperature difference between the ends of a conductor with conductanceCij.One end of the conductor is nodeiwith temperatureTiwhile the other end is nodejwith temperatureTj.Therefore,one of the keys to computeλeffusing TCNM is to calculate theCijbetween nodesiandj.

3.2.Computation of thermal conductance

Green,purple and blue cylinders in Fig.5e were used to calculate the thermal conductance in particles,interparticle contacts and near-contacts.The thermal conductanceCof a cylinder could be simply computed using the following equation:

whereλis the thermal conductivity,Ais the cross-section area of the cylinder,andLis the length of the cylinder.Then,the main task of computing the thermal conductance for different cylinders was to achieve the associatedλ,AandL.The following expressions were used to compute local thermal conductance:

(1)Particle conductance(green cylinders in Fig.5e):

whereλsis the thermal conductivity of the solid,LPis the distance between particle centroid and the corresponding interparticle contact or near-contact,andApis the cross-section area of particle cylinder and computed using a fraction coef ficientχto consider part of the particle volumeVPto avoid the overestimation of the heat transfer inside the particle.For sphere packings,0.25 was assigned toχ(Yun and Evans,2010).

(2)Thermal conductance between two particles with a gap(blue cylinders between particlesiandjin Fig.5e),which is the sum of thermal conductanceof a group of gap cylinders:

whereλvis the thermal conductivity of the void phase,Lvis the pixel size of XCT images,andis the length of each blue nearcontact cylinder in Fig.5e and computed using the distance between boundary voxels in neighbouring particles.The maximum value oflisε,which is also the cut-off to determine the effective zone for particle-fluid-particle heat transfer.

(3)Thermal conductance at interparticle contact when two particles share the same boundary voxels(purple and green cylinders between particlesjandkin Fig.5e):

whereACis the interparticle contact area computed as the sum of the areaof the shared boundary voxels.To overcome the aforementioned partial volume effect,was corrected using Eq.(6)by applying penalty parameterτto each local greyscalegnof the shared voxels and their maximum greyscaleκis another penalty coef ficient to reduce the overestimation because of neglecting roughness at interparticle contacts.It was set as 0.75 for the original greyscale XCT images since the work of Askari et al.(2015)showed that a 25%overestimation could occur if the images were left untreated after threshold segmentation.LCis the length of the interparticle contact cylinder,which is assumed to be 3Lv(van der Linden et al.,2021)referring to the studies of Bauer and Schlunder(1978)and Shapiro et al.(2004):

Once the local thermal conductance was available,the effective thermal conductance between two particles was calculated as follows for the contact condition that a gap smaller thanεexists between two neighbouring particles:

For the interparticle condition,the effective thermal conductance is computed as follows:

3.3.Computation of effective thermal conductivity

According to the equations in the previous section,the input parameters in TCNMare summarised in Table 2.For mineral grains in soil,their thermal conductivity is usually set as 3 W/(m K)while 0.025 W/(m K)for air(Young et al.,1996).The determination ofχ andεis based on:(1)the calibration of theλeffof same sphere packings from our TCNM with that from sphere packing thermal network model developed by Yun and Evans(2010);and(2)the best match with theoretical and experimental results shown in Section 5.1.A quarter of particle volume was used to compute the λeffof sphere packings,i.e.χwas 0.25 in Eq.(3)according to Yun and Evans(2010),and in this work,a value of 0.23 was assigned toχto best match the theoretical and experimentalλeffof natural sands made of irregular particles.The value of the cut-off rangeεof near-contact in a sphere packing was recommended as the harmonic mean,i.eof two neighbouring particle radiiriandrj(Yun and Evans,2010),we used the average,i.e.0.5(ri+rj),for natural sand after calibration processes(van der Linden et al.,2021).The original greyscale XCT images were used to compute theλeffof original sands,andκ=0.75 andτ=10 were imported to TCNM to reduce the overestimation during the transmission from greyscale to binary images(Fig.5a and b)(van der Linden et al.,2021).As binary images were used to compute theλeffof cemented samples,penalty was not required forκandτto deal with the binary images.Hence,they were set as 1 and 0 for the computation ofλefffor cemented samples,respectively.

Table 2Input parameters in the thermal conductance network model.

After setting the required input parameters,boundary temperatureTtop=293 K andTbot=292 K were imposed to the top and bottom nodes in the thermal network,respectively,as shown in Fig.5f,to generate a thermal gradient.Next,a Python library OpenPNM(Gostick et al.,2016)was used to compute the temperature and heat flux at each node by importing the local conductance and boundary temperatures to Eq.(1)followed by successive iteration.Finally,the effective thermal conductivity of the sample can be computed as

where∑Qijis the sum of the heat flux of all the nodes on the crosssection.In this paper,both the top and bottom planes were selectedto compute their ownλeffand their average was treated as theλeffof the sample.

4.Results and discussion

4.1.Model veri fication

For the veri fication of the TCNM,theλeffvalues of glass beads and angular sand at rest condition as well as Ottawa sand under four stages of uniaxial compression(0 MPa,2 MPa,6.1 MPa,and 10.2 MPa)were computed using both TCNMand FEMimplemented in COMSOL Multiphysics(COMSOL AB,2020).The steps of calculatingλeffof natural sands based on XCT images using the FEM is detailed elsewhere(Fei et al.,2019b).Theλeffvalues of Ottawa sand and angular sand at rest were measured by the authors using a thermal needle probe(KD2 Pro thermal properties analyser from Decagon Devices,Inc.)following the standard ASTM D5334-14(2014).Theλeffresults measured by Yun and Santamarina(2008),Narsilio et al.(2010b),and Tarnawski et al.(2011)and calculated using the theoretical formulations are summarised in Table 3,and comparison was also conducted,as shown in Fig.6.

Table 3Theoretical models for predicting effective thermal conductivity.

To show the difference between the theoretical formulations,the variation ofλeffwith porosity from 0 to 1 is shown in a semi-log plot in Fig.6a.λeffis the same as the solid thermal conductivity when the porosity is 0,and same as the void thermal conductivity when the porosity is 1.Results from FEM,TCNM and experiments sit within the theoretical results.Speci fically,FEM results are located between the self-consistent(SC)and geometric mean(GM)models,which are larger than the TCNM and experimental results located between GM and the lower bound of the Hashin-Shtrikaman(HS-L)model.λeffobtained from TCNM shows a good agreement with experimental results in Fig.6b.In contrast,FEM results are larger than experimental results,presenting obvious overestimations,especially when the porosity is small.The decrease of porosity means more contacts and stronger possibility of overestimating interparticle contact area,leading to the overestimation ofλeff.Similarly,in the uniaxial compression tests of Ottawa sands,the overestimation increases when the sample is subjected to a larger loading.

4.2.Effective thermal conductivity of surface cemented sands

Apart from applying an external loading to enhance theλeffof original sands,as shown in Fig.6,the interparticle contact area can be directly increased by cementation to obtain a largerλeff.For each sand,the original uncemented sample and 11 cemented samples were prepared and theirλeffvalues were computed using TCNM.For Ottawa sand,itsλeffincreases from 0.29 W/(m K)to 0.44 W/(mK)after applying 10.2 MPa axial stress,as shown in Fig.7(particle or cementation breakages are not accounted for).In contrast,extra solid caused by cementation only occupies the gaps between real interparticle contacts and makes the point-to-point contacts become surface-to-surface contacts in D0 sample(i.e.partially filled void voxels become solid phase in Fig.3).The increased contact area enhances theλeffof the original Ottawa sand from 0.288 W/(m K)to 0.588 W/(m K).In D2 sample,a thin surface cementation was generated,as shown in Fig.2d,which increases theλeff150%.One can achieve a reduction in porosity(compaction or densi fication)in a site by the traditional preloading ground improvement technique,but since 10 MPa stress is not easy to be achieved using traditional preloading on an engineering site,it follows that cementation could be a more accessible approach to enhance theλeffof soils in practice(and would not prompt particle breakage).

Fig.6.Comparison of effective thermal conductivity among different resources from FEM,TCNM,experiments and theoretical models:(a)The theoretical results in log scale with porosity ranging from 0 to 1;and(b)A subplot of(a)with a smaller range of porosity andλef f in normal scale.Note that the colours of markers present different sands while the shapes of markers indicate how theλef f was obtained.

Fig.7.Comparison of the effects of cementation and compression onλef f from the thermal conductance network model.

In the original sands,angular sand has the largest porosity and the least solid among the three sands,as shown in Fig.8.However,the porosity decreases and extra solid volume increases at the fastest speed in angular sand,which later leads to a largerλeffthan that of glass beads and Ottawa sand at the same porosity.

The cemented samples D0-D10 for each sand are related to the cementation level generated by pixel dilation.The correlation between cementation level and the ratio of theλeffof each cemented sand to that of the original sand is presented in Fig.9a.Angular sand has the highest increase inλeffat each cementation level,which indicates that theλeffof a sand consisting of a more complex particle shape may be more prone to be thermally enhanced.The faster increment ofλeffin angular sand than that in others is owing to the quicker increase of its average coordination number(CNave),as shown in Fig.9b.The coordination number corresponds to the interparticle contact rendered as red edges in Fig.5e,acting as an essential heat transfer path in granular materials(Yun and Santamarina,2008).

4.3.Effective thermal conductivity of directional cementation

Instead of injecting cement into the soil to reinforce the ground,MICP is another engineering method to enhance cementation between sand particles(Wang et al.,2020),and the cementation may be anisotropic due to the fluid flow direction and pressure or temperature gradient.The anisotropic cementation can intensify the thermal anisotropy of the sands since the internal preferential heat transfer path could be changed.

For both vertically and horizontally cemented samples,heat transfers through themvertically,as shown in Fig.4.Hence,vertical cementation enhances the particle connectivity along the heat transfer direction,which enforces the heat transfer in a shorter path,more straight from top to bottom.Hence,the vertical cementation results in a higherλeffthan that in the horizontal direction,as shown in Fig.10,even though at the same porosity,which also indicates the importance of microstructure for heat transfer more than mere porosity.The difference betweenλeffdue to cementation direction also offers a guide to the engineering application:calcite is desired to be generated along the heat transfer direction when using MICP to improve theλeffof the ground.As for the three sands either after vertical or horizontal cementation,theλeffof angular sand is still improved to the most extend,similar to the surface cementation in Section 4.2.Hence,improving λeffin irregular sands by cementation is more signi ficant.

Fig.8.λef f from the thermal conductance network model versus porosity in cemented sands.Note that the x-axis is inverted to present the increase of solid volume due to cementation.

Fig.9.The effect of cementation level on the increase of normalisedλef f from(a)thermal conductance network model and(b)normalised CN.

5.Conclusions

In order to improve the ef ficiencies of geothermal systems,surface and anisotropic cementations were imposed to three sands varying in particle shape,aiming to achieve superior effective thermal conductivity(λeff).A TCNM was implemented to calculate theλeffof the sands.The comparison betweenλeffof uncemented sands from TCNM,FEM simulation and experiments demonstrates thatλefffrom TCNM aligns with experimental measurement and TCNM can overcome the potential overestimation ofλeffin FEM simulations.The cementation of the sands in this paper was simulated by increasing the mineral fraction gradually,thus the interparticle contact area increases at the same moment.Theλeffof the cemented samples from TCNM fall within theoretical results,showing the robustness of TCNMfor simulating granular materials with large contact areas.

Cementation can increase the interparticle contact area more ef ficiently than applying external pre-loading of the ground.Additionally,a sand with more irregular particles show a larger enhancement ofλeffduring cementation than more rounded particle sands.It is also noticeable that cementation along the heat transfer direction can create more preferential heat transfer paths and thus lead to higherλeffthan that of the condition when cementation is perpendicular to heat transfer.

Fig.10.λeff of three sands after directional cementation.

Declarationofcompetinginterest

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to in fluence the work reported in this paper.

Acknowledgments

Project ARC DP210100433,and Dr.Anton Maksimenko and other beam scientists of the Imaging and Medical Beam Line(IMBL)at the Australian Synchrotron are acknowledged for their support.