Thermal Cracking Analysis of Microbial Cemented Sand under Various Strains Based on the DEM

Microbial-induced calcite precipitation (MICP) is a novel ground improvement method to increase the strength and stiffness of sand. However, the influences of temperature load on the internal microstructure of microbial cemented sand (MCS) material under the experimented strain have always been a key concern for the extensive application. *ree kinds of experiments, X-ray diffraction (XRD), X-ray computed tomography (XCT), and scanning electron microscopy (SEM), were conducted to explore the composition, shape, and bonding characteristics of physical assemblies in this paper. A precision DEMmodelling of MCS, mainly composed of irregular particle modelling and a mesoparameter calibration algorithm, has been proposed for the thermal cracking analysis under various strains (i.e., 1.0‰–3.0‰). Research results indicate that three kinds of bonding (that is sand-calcite, calcitecalcite, and sand-sand) are present in the MCS material. *e application of temperature has a superposition effect on the damage of MCS material with increasing strain. Moreover, as the heating duration gradually increases, the effect of thermal rupture produces a distinct quiet period.*e length of thermal cracks in the transverse direction increases throughout the heating process.


Introduction
Microbial-induced calcite precipitation (MICP) is a novel ground improvement method that is used to increase the strength and stiffness of sand using natural biogeochemical processes. Currently, MICP technology is being gradually applied in the field of engineering [1][2][3]. As a reinforcement technology, microbial cemented sand (MCS) structure based on MICP technology will inevitably encounter an extremely harsh geological environment, such as high stress and high temperature, which will result in the failure surface at the internal yield point of the internal structure as the strain increases [4]. If the deformed microbial cement structure is exposed to a temperature load, the internal structure of the MCS material might be affected. At the mesolevel, various internal microdefects may occur and develop; finally, it will change the basic mechanical properties of the MCS structure. It is difficult to observe the development of the microcrack network of the structure under the influence of temperature from the macroscopic physical experiments. e discrete element method (DEM) models are wellsuited for providing the necessary data to analyse the interactions in granular mixtures and are beneficial for improving our understanding of thermal cracking behaviour [5]. Recently, many studies have used the DEM method to perform thermal cracking analysis of rock-like materials [6][7][8][9][10][11][12]. However, the analysis of the fracture mechanism of deformed cemented structures under temperature load and, in particular, MCS material is still lacking, which is mainly due to multiple interactions along the interface among different particles (i.e., sand and calcite) in the MCS material [13]. is behaviour leads to the complexity of establishing both mesoscale bonding model and coupled thermomechanical model.
is study is aimed at providing an improved understanding of the thermal cracking behaviour of MCS material on a microscopic scale, which is not easily elucidated by macroscopically experiments. e remainder of the following sections has been organized as follows: First, the physical experimental process of MCS has been introduced briefly. e results obtained using X-ray diffraction (XRD), X-ray computed tomography (XCT), and scanning electron microscopy (SEM) were analysed to explore the material physical characterizations at different scales. Second, the observed evidence was used to set up the mesoparameters of the MCS material model with the aid of the adaptive differential evolution (ADE) algorithm; the interactions of sand-sand, sand-calcite, and calcite-calcite were carefully considered. Based on this work, the extensive coupled thermomechanical model was finally constructed to model the thermal cracking behaviour of the MCS material, with five different strains (i.e., 1.0‰, 1.5‰, 2.0‰, 2.5‰, and 3.0‰) and various heating durations (i.e., 0.02 s, 0.05 s, 0.10 s, 0.25 s, and 0.50 s) being considered.

Experimental Setup and Experimental
Characterization Analysis of Physical Assemblies

Experimental Setup.
e experimental MICP treatment process and uniaxial compressive strength and mineralogical and textural analysis tests of the MCS material are summarized below.

Sample Preparation.
e basic physical properties of standard sand used in this study are shown in Table 1. e paper briefly introduces the sample preparation method; however, the detailed method has been described by Xu et al. [14].
A visible light spectrophotometer (WFJ2000, UNICO) was used to detect the number of microorganisms. e wavelength was 600 nm, and the measured value was OD 600 . A cementation solution of urea-CaCl 2 was used where the urea was the nitrogen and energy source of microbial growth and CaCl 2 was the calcium source in the process of MICP. e test conditions selected in this paper were OD 600 � 1.50 and [Ca 2+ ] � 0.75 mol/l. Five additions of cementing solutions and one of the bacterial liquids were considered one cycle for the entire cementation process, and two cycles were performed for the test sample.

Uniaxial Compressive Strength (UCS) Test.
e UCS testing of the MCS material was conducted using a microcomputer-controlled electrohydraulic servo loading testing machine (TAW-2000) (shown in Figure 1). e maximum loading capacity of TAW-2000 was 2000 kN. e uniaxial compressive strength test was set to be displacement controlled at a loading rate of 0.05 mm/min until failure occurred. A Micro-II digital acoustic emission (AE) system was used simultaneously to identify fracture and deformation events caused by damage to the sample.

Mineralogical and Textural Analysis Test.
e mineralogical and textural analysis test consisted of three parts: XRD, XCT, and SEM. e cemented sand body on the surface of the column was examined using a scanning electron microscopy (SEM) and X-ray diffraction (XRD). e paper gives a brief introduction to the methods; however, detailed methods and relative parameters have been described by Xu et al. [14].
A micro-XCT instrument (diondo d5) with the resolution of 27 µm, a voltage of 180 kV, and tube current of 0.18 mA was used, and the exposure time of the detector was 1000 ms in this study. e scanning process (shown in Figure 2) consisted of four parts: (a) test sample selection, (b) placing and adjusting the sample position, (c) setting the parameters and scanning the sample, and (d) feature extraction.

Composition of Physical
Assemblies. e XRD patterns of the pure sand and cemented sand are shown in Figure 3. According to the XRD patterns, the pure sand was composed of quartz, and the mineralized cemented sand was composed of quartz and calcite phases.

Shapes of Physical Assemblies.
e XCT experimental result is shown in Figure 4. Obviously, the calcite crystals in the MCS samples were irregularly shaped and randomly distributed, and most of the particles were elongated and polygonal.

Bonding Characteristics of Physical Assemblies.
According to the results of the microscopic images shown in Figure 5, it was found that three kinds of bonding were inside the cemented sand: sand-sand bonding, sand-calcite bonding, and calcite-calcite bonding. e sand-sand bonding exhibited an inhomogeneous bonding between particles ( Figure 5(b2)), with some visible holes at the interfaces of the bonding, while a tight bonding was exhibited for the sand-calcite particles and calcite-calcite particles ( Figure 5(b1)). Considering the aforementioned bonding characteristics of the MCS material, it is necessary to take into account the different bonding characteristics and distribution of the MCS material when setting the mesoparameters of the numerical model.

Numerical Specimen Generation and Parameter Calibration
In the current study, although the PFC 3D version could be employed for solving 3D problems, the PFC 3D simulations required very long computational times, and the boundary conditions were difficult to control [15]. On the contrary, the 2D simulations have the advantage of explaining key phenomena and mechanisms in less time than 3D modelling. In addition, the visual deformation patterns and force chains are easily captured by 2D simulations. Many previous 2 Advances in Materials Science and Engineering  studies [16][17][18][19] illustrated that the 2D simulations can adequately capture various complex mechanical features of the materials when considering the shape's e ects. erefore, the numerical simulation for biocemented sand specimens was carried out using a PFC 2D model in this research.

Numerical Specimen Generation.
e numerical specimen generation process was conducted through the following steps.

Contact Constitutive Model Selection.
According to the arrangement of MICP in the microscopic images in Figure 5, the bond pattern was close to the parallel-bond model ( Figure 6) of DEM [20].
e parallel-bond model (PBM) was depicted as a rectangle of cement-like material, which could transmit both force and the moment between the particles [21]. Moreover, the PBM has been successfully used to explore the mechanical behaviour of MCS material [22]. As such, the PBM was selected for the numerical study in this research.

Calcite and Sand Particles
Modelling. As in the aforementioned XCT scanning results, the calcite crystals in the MCS samples were irregularly shaped and randomly    distributed, and most of the particles were elongated and polygonal. ere were a larger number of physical contacts for irregular particles than that for circular particles [23]. Moreover, the particle shape played a key role in the behaviour of the material [24]. To overcome these problems, the appropriate geometrical shape of irregular particles should be modelled. is study adopted an irregular generation algorithm to determine the shape e ect of the MICP material. e irregular polygon generation process of calcite particle is shown in Figure 7, which was described in detail by Tang et al. [13].
For the sand particles, the previous analysis results [14] have shown that the particles were mostly circular, and the particles for this particle size have been characterized by circular particles [25]. For this reason, the sand particles were simulated using circles.

Numerical Specimen Generation.
Using the parallelbond model, a 2D numerical model with the height of 100 mm and the width of 50 mm was rst built to simulate the laboratory behaviour of a cylindrical specimen of biocemented sand with the same size as that of the numerical model. According to the research results [13], the calcite content of cementation sand was 3% in the 2D model, and the numerical model consisted of 3% calcite and 97% quartz by area, combined with the generation method of calcite.
Subsequently, the numerical particles were produced in the numerical specimen model based on the grain-size distribution (GSD) of the real particles. e numbers of quartz and calcite particles were 165244 and 33955, respectively. e GSD for real particles and numerical particles are shown in Figure 8; a strongly consistent GSD was exhibited for real particles and numerical particles.

Parameter Calibration.
Realistic mesoparameters were required when using DEM to calculate the numerical value. However, the mesoparameters of the numerical model were generally adjusted by trial and error until the numerically simulated results corresponded to the macroparameters of the real specimens.
is method usually requires a large number of trial calculations, and the calculation period is long; thus, the human error can be signi cant. Compared to the trial-and-error method, adaptive di erential evolution algorithms (ADEs) have been widely used in many elds and have the advantages of simple control parameters, global search ability, strong robustness, fast convergence speed, and precision characteristics [26][27][28]. ey would also be suitable for solving the problem of parameter inversion, especially for the inversion of discontinuous, nonlinear, nondi erentiable, and analytic functions.

Calibration Method.
In this paper, the mesoparameters of MCS material were calibrated using the ADE method. e calibrated mesoparameters include E * , E * , k * , k * , fric, pb_ten, pb_coh, and R.
e SEM characterisation of the material (shown in Figure 5) makes it necessary to consider the di erent bonding characteristics and discreteness of MCS material when setting mesoparameters. According to the results of Feng et al. [22], di erent bonding characteristics can be characterized by di erent bonding radii. When the parallel bonding radius, R (m), is determined, the particle ow theory suggested that the mesoscale bonding radius, R, was related to the radial coe cient, λ. For the sand particlecalcite and calcite-calcite bonds, λ was assigned a xed value of 0.5 because these bonds were relatively strong according to the SEM analysis results. For the sand particle-sand particle bond, a random number to vary the parallel-bond radial coe cient λ between 0.25 and 0.5 was introduced based on the previous PFC study on microbial cemented sand column by Feng et al. [22]. us, the bonding radius was determined and is not included in the parameter inversion procedure.
e mesoparameters M of the MCS material used in the simulation can be expressed as where E * is the e ective modulus (force/area), E * is the e ective modulus of the bond (stress), k * is the normal-toshear sti ness ratio, k * is the normal-to-shear sti ness ratio of the bond, fric is the friction coe cient, pb_ten is the tensile strength (stress), and pb_coh is the cohesion (stress). e adaptive di erential evolution algorithm used to invert the mesoparameters is shown in Figure 9, and the overall process is described in the following three steps.
(1) Construction of the Characteristic Parameters. Using the PFC simulating uniaxial compression test, N i (i ≥ 20) groups of mesoparameters M j [j ∈ (1, N i )] were constructed, and the stress-strain and crack data corresponding to these parameters were obtained. It was noted from the laboratory test results of elastic modulus size of the MCS material that the elastic modulus of sand particles is relatively large and that of calcite is relatively small [29].
us, the initial elastic modulus parameters were set according to this size. Subsequently, based on the stress-strain and crack data, the characteristic parameters were calculated using MAT-LAB (R2015b) software. e obtained characteristic parameters were used as the initial population samples

Advances in Materials Science and Engineering
, elastic modulus (E j ), and the crack initiation stress (σ n−j ). Note that the crack initiation stress (σ n−j ) was identi ed as the axial stress where the number of cracks reached 1% of the total number of cracks at peak stress [20]: (2) Determination of the Target Inversion Value from the Uniaxial Compression Data. e experimental uniaxial compression tests provided the real stress-strain behaviour and AE monitoring data. According to the stressstrain relationship and the de nition of the cracking stress, the peak stress (σ c ), peak strain (ε c ), elastic modulus (E), and crack initiation stress (σ 1−p ) of the MCS material were obtained and used as the target values for inversion. It should be noted that AE monitors lowenergy high-frequency seismic events resulting from crack initiation, so it can be used for the identi cation of σ 1−p [30]:    Advances in Materials Science and Engineering (

3) Construction of the Fitness Function of the Inversion
Algorithm. According to the combination of di erent mesoparameters M j [j ∈ (1, N i )] , the study focused on the relatively small errors between the simulated and experimental characteristic parameters where F is the target value. e objective function of the inversion of the mesoparameters can be de ned as Furthermore, the tness function was constructed as follows: It was found that the larger individual tness values resulted in smaller objective function values; the mesoparameters re ected the real MCS material properties.

Calibration Results and Comparison of Calibration
Methods. According to the above calibration method, the mesoparameters with small errors were observed, as shown in Table 2.
e uniaxial stress-strain curves and failure characteristics calculated from the numerical simulation are compared to that of the experimental data in Figure 10.
e results shown in Figure 10 were found by comparing the numerical results to the results of physical experiments, with 1.95% and 4.0% as the relative errors in the peak stress and peak strain, respectively. Additionally, the numerical tests revealed a similar trend between the macroscale cracking and the overall fracturing feature with minor local deviations. e deviations are primarily due to the discrete connection of the microscale cracks in the macroscale cracks, which are inherently induced by the discrete contacts between particles in the model developed in this study.
Meanwhile, some simulation results di er from the testing results; speci cally, the elastic module obtained in each case di ers from the others. is di erence is justi ed by the model reconstruction of microscale particle structure and the related physical properties. In a numerical simulation, it is almost impossible to reconstruct a model sharing the same mesoscale structure feature and physical properties as those of a physical testing sample [31]. In view of these factors, the mesoscale parameters used in this study are appropriate for characterising the MCS material properties.
Currently, the trial-and-error method has been the most commonly used calibration method for selecting appropriate mesoscale parameters [32,33]. e stress-strain curve corresponding to the mesoparameters obtained using the trial-and-error method was plotted on the same graph ( Figure 10) in order to compare the results of ADE and the trial-and-error method. Figure 10 clearly illustrates that the numerical curve obtained by the ADE method is similar to the curve obtained by Tang et al. [13] using the trial-anderror method. By comparison with the mesoparameters obtained by Tang et al. [13], most mesoparameters are almost the same, except for pb_ten and pb_coh. e mesoparameters used for the trial-and-error method usually have to be assigned according to the mesoscopic properties of the material. For relatively simple contact types, the trial-and-error time may be relatively short. However, when the material involves multiple contact types, the trial-and-error process is extremely complicated and takes more time. Compared to the trial-and-error method, ADE is able to create an initial database and obtains target parameter values through parameter inversion. e entire process only takes a certain amount of time during the initial database setup process. Once the initial database is established, corresponding mesoscopic parameter values can be quickly obtained. Furthermore, the established database can provide initial data for the subsequent acquisition of mesoparameters of related materials with similar mechanical behaviour.

Coupled Thermomechanical Model
In the following section, the thermal-mechanical model is constructed.

ermal Model. e heat conduction in a continuum is given by
where q i is the heat ux vector (W/m 2 ), q v is the volumetric heat intensity (W/m 3 ), ρ is the density of the material (kg/m 3 ), C is the heat capacity under constant pressure (J/kg·K), and T is the temperature (K). e heat out ow per unit volume is given by the divergence of q i . e average value of the divergence of q i in the reservoir is given by e volume integral can be replaced with a surface integral by the application of the Gauss divergence theorem to the reservoir. In addition, if the heat ow is restricted only to the N thermal pipes associated with the reservoir, the surface integral can be replaced by the following summation: Hence, If each pipe p is regarded as a one-dimensional object with a thermal resistance per unit length of η, then the power in a pipe would be given by where ΔT is the temperature di erence between the two reservoirs on each end of the pipe and l (p) is the pipe length. A heat source, Q (p) , and thermal resistance, η, are associated with each pipe, which has a nite length equal to the distance between the centroids of the two particles that it connects. For a two-dimensional model, the thermal resistance, η, is given by where n is the porosity of the numerical specimen, V (b) is the volume of the ball (b), l (p) is the length of the thermal pipe p, N b is the number of balls, K is the thermal conductivity, and N is the number of thermal pipes. For a two-dimensional system, the ball disk in the PFC is unit 1. Hence, V (b) S (b) . From equations (6)-(11), the governing equation can be rewritten as follows: e amount of thermal energy deposited into a material due to microwave heating, which is equated to the volumetric heat intensity, is given by [7] where P d is the power density (W/m 3 ), f is the frequency of the microwave irradiation (Hz), ε 0 is the permittivity of free space (8.854 × 10 −12 F/m), ε r ″ is the dielectric loss factor of the mineral, and E 0 is the electric eld intensity resulting from the microwave power (V/m). Substituting P d in place of Q v in equation (12) gives us, the temperature increase can be obtained by solving equation (14) for a given power density and heating duration.

ermal-Mechanical Coupling.
e temperature field can be simulated using the governing equations for microwave heating. Based on the definition of different material properties, the conversion of microwave energy was achieved by setting the power density and heating duration in PFC modelling.
When the microwave energy was converted into heat energy, thermal strains were produced in the material, resulting in changes in the particle radii and the force transmitted by each parallel bond [34]. is should be taken into consideration in the heating of both the particles and the bonding material. For a temperature change ΔT, the change in the particle radius R is given by where α is the coefficient of thermal expansion of the particle.
If a parallel bond is present at the contact associated with a pipe, then where K n is the normal bond stiffness, A is the crosssectional area of the bond, α is the expansion coefficient, L is the bond length (assumed to be equal to the distance between the centroids of the two particles connected by the pipe), and ΔT is the temperature change.

Boundary Conditions and ermal
Parameters. e definition of the boundary conditions generally requires a description of the local boundary in terms of an outward normal unit vector. More realistic boundary conditions also consider the convection and radiative boundary properties between the specimen and the surrounding air [35]. However, a material modelled in the PFC 2D software has no explicit representation of its boundary (unlike in a continuum code, where zone surfaces comprise the boundary). is problem is overcome by specifying the power input to a boundary surface instead of the flux. e external boundaries of the thermal model were considered to be thermally insulated. erefore, the convection and radiation were not incorporated into the model since the microwave heating durations were quite short.
Considering that the thermal input properties of an MCS was determined by the microstructural components, two different sets of thermal parameters were obtained from the literature [36][37][38], as presented in Table 3. Based on a previous study [33], the microwave heating durations were between 0.02 s and 0.50 s, and the power density was equal to 1 × 10 9 W/m 3 in this study.

Heating Effect Analysis Using Coupled
Thermomechanical Model

Heating Effects.
e temperature fields as a function of the MCS material under different heating durations are shown in Figure 11, where a gradual increase in the sample temperature as a result of longer simulation time of heating can be observed. e distributions of different particles at the computational time of 0.50 s in the testing sample were illustrated with sand in crimson and calcite in light red. It should be noted that the heating-induced cracks in the sample have not been accounted for in this study. e focus is on the influence of short-term heating on the deformed samples under different magnitudes of strains.

Influence of Strain Magnitudes on the Failure
Characteristics.
e results of the simulation using the thermomechanical model were interpolated using a python language compiler with the uncoupled model being included for a direct comparison ( Figure 12). is article defines two types of cracks: cracks that occur only under mechanical loading are called mechanical cracks (shown in blue), and cracks that occur under temperature loading are called thermal cracks. It should be noted that the same heating duration of 0.50 s (corresponding to T � 210°C) is adopted under five different strain magnitudes.
As shown from Figure 12(a), the results captured from uncoupled modelling can be summarized as follows: (1) when the strain was 1.5‰, few discrete microcracks were present throughout the sample; (2) when the strain exceeded 1.5‰, there were mechanical cracks in the structure (inside there exist damage) and mainly in the quartz matrix; (3) under larger strain (3.0‰), more microcracks gradually commonly coalesced, with few cracks in the calcite.
As shown from Figure 12(b), thermal influence on the cracks can be observed; when the strain exceeded 1.5‰, the number of microcracks (shown in cyan) caused by thermal cracking continuously increased. Specifically, it was found that a large number of microthermal cracks (shown in the green frame) formed and coalesced in the mechanical crack tip (shown in blue) and its surroundings when the strain was 3.0‰. ese indicated that the mechanical crack (shown in blue) acted as a "rock bridge" with the increasing strain, where the damage to the material played a guiding role in the thermal crack (shown in cyan).
Besides, as shown from Figure 12(b), the thermal crack in the material mainly extended to the quartz matrix. Although the hardness of quartz is clearly higher than that of calcite, the bonding between the quartz particles in the MCS material is determined by the cementation. Compared to the calcite particle ( Figure 5), weak cementation occurred in the quartz sample (some holes were observed at the bond interface that decreased the cementation strength and degraded the cementation performance), resulting in uneven bond strengths. erefore, the bond rupture was more likely to occur, and the crack initially propagated through the quartz matrix.
According to the simulation results shown in Figure 12(b), we can also see that the crack surface was rough during crack propagation and that the crack bifurcation was easily observed when the strain was 3.0‰. For the crack bifurcation phenomenon, it showed that, under the action of stress, the crack propagation path with the minimum impedance is taken [39]. is demonstrates the heterogeneity of the material and confirms that the numerical model is reasonable for considering the strength distribution of the actual sample.

Influence of Increase in Temperature on the Failure
Characteristics. Figure 13 shows that thermal cracks (shown in cyan) for a sample, with a strain of 3.0‰, were observed as the temperature increased. During the heating process, a tensile stress caused by microwave energy is formed in the material [40][41][42]. At a relatively lower temperature, the tensile stress generated is smaller. e fracture effect for this thermal crack developed a quiet period with an increase in time for the duration of 0.25 s (corresponding to T � 110°C) of the quiet period. With the increase in temperature, the tensile stress rapidly increased after the experience of higher temperature, and thermal cracks gradually converged and interacted, resulting in the appearance of obvious failure zones (Figure 13: 0.50 s (corresponding to T � 210°C)).
ese results were obtained because the increase in heating duration gradually raised the temperature of the material. e temperature in the MCS material was not very high for a shorter heating period, and the tensile stress generated potentially promoted the fracture at a location of the model with lower internal bond strength. A location with higher bond strength was not sufficient to induce bond fracture. erefore, in the short term, the fracture of the model was not very obvious. Subsequently, the temperature rose rapidly as the heating duration increased, and the generated tensile stress exceeded the bond strength between the particles. e new thermal cracks developed rapidly, and many small thermal cracks gradually converged and interacted, resulting in the appearance of obvious failure zones (shown in the cyan lines in the final sample). Finally, the thermal cracks exhibited the behavioural characteristics across the structural plane in the transverse direction.
In summary, the temperature gradually changed with an increase in heating duration during the heating process, and the material inevitably produced a rupture, which further led to a change in the macromechanical behaviour of the MCS material. e application of temperature has a superposition effect on the damage of the MCS material.

Comparison between Calculated and Experiment Results.
As mentioned above, for the MCS material, the cracks were preferentially located in the internal damage area of the structure under the action of temperature. With the increase of temperature, many small thermal cracks gradually converged and interacted with the obvious failure zone. To validate the accuracy of numerical simulation results, a heating test was carried out. e heating test has been simplified to verify the resulting effect of temperature on the existing damage zone. e experimental specimens were selected after uniaxial compression, and the heating experiment was conducted using heating equipment (KW-202). e longest heating duration, 285 min, was chosen to allow the temperature to reach at approximately 210°C, which corresponded to the maximum value in the numerical analysis. Subsequently, the damage images of the MCS material under the different heating durations via 3D Super Depth Digital Microscope (VHX-5000) were observed. Figure 14 shows the heating test process with the test results shown in Figure 15.
From Figure 15, it can be observed that when the heating duration was 15 min, the particles in the damaged area were closely bonded, with a few local cracks (Figure 15(b)). As the heating durations increased, more localized cracks were observed, and the interaction of the localized cracks occurred (Figure 15(e)). When the heating duration reached 285 min, the crack had spread from the damage boundary to the surrounding area (as shown in Figure 15(f )). According to the test results, the failure mode of the experimental results is consistent with the simulation results, and it confirms that the simulation results are reasonable.

Conclusions
e material physical characterizations at different scales were analysed in this paper. A method for calibrating the input MCS material properties was designed, and a coupled thermomechanical model was constructed. e behaviour of the MCS material subjected to temperature loading was investigated under strains varying from 1.0‰ to 3.0‰. e following conclusions were drawn from the study: (1) From a comparative analysis between ADE and trialand-error method, the adaptive differential evolution algorithm was successfully used to explore the MCS mesoparameters and can also be used to establish the comprehensive database of parameters for the MCS material in DEM analysis. (2) For the wide range of the applied strain magnitudes (lower than that corresponds to the failure state), the elastic deformation can be observed under low-strain magnitudes (<1.5‰), while the progressive development of the localized cracks is more evident under high-strain magnitudes. ese cracks occur more frequently at the sand-sand interface due to lower bonding strength among them.
(3) For the typical strain of 3.0‰ associated with the prefailure state under the different temperatures, thermal cracks behaviour experience two phases. (1) A quiet period: local crack progressively developed and the thermal effect on the crack might be ignored. (2) An explosion period: the local cracks gradually converged and interacted within the distinct failure zone. e length of the thermal cracks in the transverse direction increased throughout the heating process.

Data Availability
e data used to support the findings of this study are available from the corresponding author upon request.

Conflicts of Interest
e authors declare that there are no conflicts of interest regarding the publication of this paper.