Effects of sample preparation on triaxial compression tests formed with uniform spherical particles

$e$

= void ratio;

$VT$

= total sample volume;

$VS$

= solid volume;

$VV$

= volume of voids;

$H$

= measured sample height;

$D$

= measured sample diameter;

$Ws$

= sample weight;

$ρ$

= sphere’s density;

$shi$

= height of slice, $i$⁠, measured from the base of the sample;

$hi,j$

= distance between centroid of particles, $j$⁠, to the base of the sample as a function of $x,y$and $z$⁠;

$di,j$

= distance between $shi$ and $hi,j(x,y,z)$⁠;

$R$

= radius of spherical particle;

$ri,j$

= projected cross-sectional radius of spherical particle, $j$⁠, on a slice plane, $i$⁠;

$Asi$

= sum of all projected areas of spherical particle on slice $i$⁠;

$At$

= cross-sectional circular area of the triaxial mould; and

$epi$

= planar void ratio at slice $i$⁠.

Understanding the mechanistic behaviour of granular materials has drawn attention from various fields from pharmaceutical manufacturing to geosciences, engineering, mining and agriculture. Numerous studies have shown that engineering behaviours of granular materials are significantly affected by sample preparation methods under controlled laboratory conditions (Oda, 1972; Arthur and Menzies, 1972; Ladd, 1974; Mulilis et al., 1977; Nemat-Nasser, 1980; Pradhan et al., 1989; Nakai, 1997; Vaid et al., 1999; Wan and Guo, 2001; Yamamuro et al., 2008; Thompson and Wong, 2008; Sadrekarimi and Olson, 2012). Sample preparation is critical in determining the response of granular materials as granular samples with similar void ratio but with different internal fabrics show different engineering response under loading.

In the case of geomaterials, Oda (1972) and Arthur and Menzies (1972), found that the stress-strain relationship of sand depends strongly on the sample preparation techniques. Ladd (1974) observed that preparation methods significantly affect the cyclic shear strength of sand. Similarly, studies of Mulilis et al. (1977) reported the effects of sample preparation techniques to the liquefaction characteristics of sand subjected to cyclic loadings. Additionally, Sadrekarimi and Olson (2012) investigated the effect of specimen preparation method on critical state behaviour of sands. Vaid et al. (1999) recommended water pluviation method over moist tamping method as samples prepared by the water pluviation method showed dilative behaviour while samples prepared by the moist tamping method showed liquefaction behaviour when subjected to triaxial compression and extension tests. Yamamuro et al. (2008) pointed out that samples prepared by pluviated methods tend to have more stable grain-to-grain contacts, which resulted in a stable behaviour, compared to dry funnel deposition and slurry deposition methods. Thomson and Wong’s study (2008) showed that both water pluviation and moist tamping methods resulted in non-homogenous samples. However, void ratio of water pluviated sample was lower only at the top quarter of the sample height while void ratio of moist tamped sample varied along the sample height.

From previous studies, factors that affect macroscopic behaviour of granular materials are preparation methods, confining stress, sample density, grain size distribution and material fabric. Common sample preparation techniques to replicate in situ conditions are moist tamping, dry funnel deposition, air pluviation and water pluviation. Different sample preparation techniques can potentially form samples with identical void ratio but dissimilar material fabrics from one sample to another and thus potentially yield different mechanical behaviours (Sadrekarimi and Olson, 2012). Variation in sample heterogeneity has a profound effect on the mechanical behaviour of granular materials. For instance, Morgenstern (2000) provided examples of geotechnical failures that underestimated the effects of in situ soils’ heterogeneity. Sun et al. (2019) concluded that particle breakage during sample preparation by layered compaction would influence the subsequent breakage during triaxial loading. Cai et al. (2018) experimentally prepared sand specimens by dry deposition and moist tamping methods and concluded that the difference in fabric anisotropy generated by the different sample preparation led to a different mechanical response under triaxial compression tests. This study adopted drop height, number of sieves and particle size, as the control variables for its experimental design.

Based on recent results from Medina-Cetina et al. (2022) and Zhu et al. (2022) where stereo digital image correlation analysis was used to trace the boundary mechanical response of sand specimens subjected to triaxial compression using the dry pluviation sample preparation technique. By setting the experimental and numerical conditions on controlled specimen uniformity, the main objective was to explore the extent of material heterogeneity and its mechanical effects when subjected to triaxial compression. This approach sets a baseline for the response of materials based on uniform cohesionless particles and the construction of uniform specimens. Further analyses to investigate strain localisation effects were also conducted by Zhu and Medina-Cetina (2022a), where 3D kinematic operators under cylindrical coordinates were applied to a subset of boundary displacement fields of triaxial sand specimen database (Medina-Cetina et al., 2022) prepared by vibratory compaction and dry pluviation. Their results demonstrate the effect of sample preparation technique on the mechanical parameters and its influence on the evolution of localisation effects in space and time, often manifested by shear, expansion or compaction bands and ultimately, the impact sample preparation has on the material’s failure mechanism.

Previously, Ng (2004) simulated three sample preparation techniques using discrete element method (DEM) for drained triaxial tests; firstly, particles were generated randomly within a cubic box then subjected to isotropic confining pressure of 100 kPa, secondly, sample prepared by dry pluviation method and thirdly, a sample prepared by tamping method with a single layer. Dai et al. (2016) investigated three means of sample preparation: gravitational deposition, multi-layer compression and isotropic compression, using the particle flow code in 2D (PFC2D) to study the DEM effect of fabric anisotropy. Yang (2018) used YADE, an open-source DEM platform, to generate cemented granular materials in 3D for a consolidated isotopically undrained triaxial compression. Hu et al. (2022) introduced a new sample preparation simulation in 2D using PFC2D to achieve sample homogeneity and using the prepared sample, managed to closely replicate experimental direct-shear results. Chang and Chen (2018) generated non-spherical particles based on the compressible packing model in 3D using particle-flow code in 3D (PFC3D) and concluded that the variation of void ratio in non-spherical samples is closer to actual sand compared to the variation of void ratio in spherical samples. Bernhardt et al. (2014) investigated the effects of sample preparation methods in the macro- and micro-scale responses in 3D DEM simulations of a direct simple shear test.

Ghosh et al. (2016) elaborated a systematic sand bed preparation by using the air pluviation technique and studied the fall height and the sample uniformity. The sample heterogeneity can affect the mechanical response of the specimen. Therefore, to validate this claim, we perform a study to verify if the sample preparation condition affects the mechanical response of the specimen and to confirm the ability of PFC3D to accurately replicate the sample homogeneity as assumed by the experimental sample preparation process.

To the best of the authors’ knowledge, there has been no preceding effort to investigate and study air pluviation sample preparation technique in both experimental and numerical setting using uniform spherical particles for a triaxial test. Recently, DEM has been used as a data generating method to calibrate both statistical and machine learning (ML) models (Liao et al., 2020). Therefore, it is imperative to rigorously evaluate the reliability of the data generation technique, as it plays a critical role in the calibration and overall performance of predictive models. Consequently, the proposed hypothesis is that samples composed of granular materials prepared following standard procedures such as air pluviation, thought to produce homogeneous structures, show in fact heterogeneous features, which may still have significant impact on its mechanistic response. To further explore the sample’s uniform condition, this work uses only uniform particles as described above.

Furthermore, existing studies have only investigated sample response from macro - perspective without considering the influence of microstructure within the specimen. The authors set out to replicate both air pluviation and radius expansion methods to study the resultant microstructure from their two distinct sample preparation techniques in PFC3D. While there have been previous efforts to study the microstructure using various imaging technique (e.g. X-ray computed tomography and scanning electron microscope) of triaxial specimens under different sample preparation techniques (Yamamuro et al., 2008; Thompson and Wong, 2008), these detailed investigations were only applied on experimental triaxial specimen and not on triaxial numerical specimens. The outcome of this study will help improving the understanding of the micro-mechanical behaviour of granular materials with respect to meso- and macro-mechanical applications. Additionally, this study will also examine the ability of PFC3D to accurately replicate the sample heterogeneity as obtained during the experimental sample preparation process setting the basis to produce data alternative laboratory experiments. The ability to rapidly produce such data is particularly beneficial for training ML models. This study adds to our experimental and numerical modelling research to leverage DEM-generated data sets to produce physics-informed neural networks such as Ortiz-Cahun et al. (2023), Duong et al. (2023), Song et al. (2023), Zhu and Medina-Cetina (2022b), Zhu et al. (2022) and Medina-Cetina et al. (2022).

With the recent surge in the popularity of ML techniques, numerous studies have sought to integrate or couple these methods with the DEM to address complex problems in granular media. Bui et al. (2019) generated porosity and sediment distribution data set, from DEM simulations of packing gravel-bed and fine sediment infiltration, to train a predictive model based on neural network. Furthermore, Liao et al. (2020) demonstrated the feasibility of using neural network to replicate the segregation behaviour of large-scale DEM simulation of blast furnace. On the other hand, Naderiet al. (2024) introduced a novel integration technique between DEM and convolutional neural network (CNN) to parallelise and accelerated computation of contact force between mapped particles. Moreover, Lu et al. (2021) also implemented CNN to accelerate DEM simulations by introducing a physics-inspired multi-scale loss function and analysing the effects of training configurations on granular flow predictions. Clearly, there is a growing body of research focused on integrating DEM with ML techniques, reflecting significant interest in this interdisciplinary approach.

This study focuses on examining the behaviour of samples formed by uniform spherical particles, which are constituted by the air pluviation method, while controlling different heights and dispersion levels to produce different specimens’ density and void ratio, to be tested under triaxial compression. Both triaxial DEMs generated using air pluviation and radial expansion methods are first simulated and then, compared to the experimental air pluviation sample preparation results. The methodology used in this study including material selection, experimental sample preparation, triaxial compression and DEM numerical model of sample preparation and void ratio analysis of both numerical and experimental samples are described in Section 2. Results of void ratio analysis, triaxial stress-strain responses and comparison of effectiveness of numerical sample preparation methods are presented in Section 3. Conclusions are presented in Section 4.

The undrained behaviour of granular materials is affected by grain size distribution, particle shape, surface friction, density and void ratio and by sample preparation methods. To concentrate solely on the effect of sample preparation method, other factors are minimised by choosing a material that has same particle size, spherical in shape and low friction coefficient. Several materials were taken into consideration for testing: chrome steel, titanium, silicon nitride, polycarbonate plastic. According to the manufacturer, Thomson (2022) and other external tests (O’Sullivan and Cui, 2009), properties of these materials are listed on Table 1 to highlight the contrast between the considered materials and guide the selection of the material. Chrome steel spheres were selected in this study due to the fact that the material is stiff with low friction coefficient and cost-effective compared with typical geomaterials. Two particle sizes, shown in Figure 1, 3 and 6 mm diameter were used in the study to obtain the effect of particle size as well.

It is hypothesised that sample preparation could condition the failure mechanism of granular samples in triaxial compression tests. The goal of the sample preparation technique is to constitute uniformly distributed cylindrical samples, so the sample is closest to the homogenous state. A previous work by Noble (2012) shows the use of air pluviation for sample preparation using steel spherical particles. The major improvement between previous work and this experiment is the inclusion of valve for flow control and the construction of PVC frame as a funnel stand. The air pluviation method reconstituted samples by drizzling of particles through a set of diffuser meshes placed on top of a fall tube connected to the specimen’s mould. Pressure was applied on the mould using a vacuum pump to keep the membrane in contact with the mould and the reading on the vacuum pump to keep the pluviated samples upright was 55 kPa. The same confining pressure was used to prepare the numerical radius expansion method discussed below. The tube had similar diameter of the mould to obtain a uniform concentration of particles’ rain across the mould. Sieves were placed on top of the tube to adjust the dispersion level of the particle flow. Steel spheres were poured into a plastic funnel whose bottom end is about 2.54 cm (1 inch) above the sieve to settle. A classic 1-inch diameter ball valve was installed at the bottom end of the funnel to control the flow. The funnel was supported by a PVC frame. Figure 2 presents the experimental set-up of the air pluviation sample preparation method described above.

The density of the sample can be increased by increasing the length of the fall tube or by changing number of sieves of size 3.35 and 8.0 mm for 3 and 6 mm diameter particles, respectively. It is predicted that the sample would achieve a denser state as the spheres were dropped from higher distance due to higher potential energy of the spheres. Also, density of samples prepared at the same drop height would decrease as the number of sieves placed on top of the tube increased because of energy loss during contacting with the sieves. The tube length varied from 30.48 cm (1 ft) to 91.44 cm (3 ft) with 30.48 cm increments. The number of sieves also changed from having no sieve to three sieves on top of the tube.

Five samples were made for each condition for repeatability purposes. Figure 3 shows a chart with the control variables considered for the experimental design. A sample prepared out of 3-mm diameter spheres is pictured in Figure 4. The dimensions of each sample (height and diameter) and sample weight were recorded to obtain the global void ratio of the sample.

A triaxial compression test based on the ASTM (2011) is performed on samples made of uniform steel spherical particles and a confining pressure of about 55 kPa was applied to the specimen using vacuum pump once the sampled settled. Figure 5 shows the triaxial test frame before compression.

The test starts as the loading frame is raised up at a rate of 1.6 mm/min, while the confining pressure is still applied to the specimen, resulting in a strain rate of 1.6 mm/min. The transducers transfer the recorded data to the computer by the use of LabView (National Instruments, 2012). Experimental results provide adequate data to determine the strength and deformation properties of a sample composed of steel particles. One triaxial compression test was performed for each sample preparation condition according to the experimental design (Figure 3).

PFC by Itasca Consulting Group and Inc (2019) was first developed by Cundall in 1971 to simulate rock-mechanics problems, then it was applied to soils problems by Cundall and Strack (1979). PFC is classified as discrete element code since it allows distinct particles that displace and rotate independently to each other and since the particles interact only at contacts or interfaces between them, while automatically recognising new contacts at the end of each calculation cycle.

A numerical model using DEM was constructed to replicate the experimental process, from sample preparation to the triaxial compression test. The goal was to achieve a target density using different sample preparation methods and then use these numerical models (i.e. specimens) to simulate a triaxial compression test under the same experimental conditions.

A numerical specimen was designed and formed using the experimental air pluviation method (Lagioia et al., 2006; Raghunandan et al., 2012). This sample preparation method is different from widely used radius expansion method in DEM modelling, where the specimen is prepared by increasing the size of the spherical particles to achieve a target porosity (Itasca Consulting Group and Inc, 2019).

The same experimental setting was reproduced by building a twin model, comprised of a funnel, a set of meshes, a tube and the specimen’s mould (Figure 6). Multiple PFC wall elements were used to form the shape of the specimen’s mould, the tube and the funnel; while the set of sieves were replicated using PFC mesh elements; and the particles were assigned a linear contact model (O’Sullivan and Cui, 2009; Itasca Consulting Group and Inc, 2019).

A linear contact model assigns a linear function of normal displacement (i.e. overlap) to normal force upon contact between particles and the normal force is restricted to be compressive only. Moreover, the allowable tangential force (i.e. shear force) obeys the Coulomb friction law and increases linearly with relative shear displacement. Slip between two adjacent particles only occur when the shear force is equal to the allowable shear force.

The experimental parameters were assigned to the PFC-3D twin model, including the wall dimensions, particle diameter, particle density, friction coefficient between particles and between wall and particle. The numerical parameters of the twin PFC-3D model, included the normal and shear modulus of particles, were chosen based on previous studies (Vaid et al., 1999; Wan and Guo, 2001; Thompson and Wong, 2008; O’Sullivan and Cui, 2009). No material properties were assigned to the PFC wall and mesh elements. Only the gravitational force was activated as an external force applied on the sample throughout the particles’ air pluviation process. Spherical particles were generated above the funnel and a wall was created to cover the discharge opening of the funnel. Then 1 m cycles were run to have the spheres deposited inside the funnel until these reach a resting condition. After that, the wall at the end of the funnel was deleted to release all the particles into the sieves, what started the particles air pluviation process. Table 2 presents a list of parameters used in the twin PFC-3D model. Pictures of the sequence described above are shown in Figure 6. The particle density, shear modulus and normal modulus were set to be identical with the physical spherical particles used in the experimental air pluviation given in Table 1. On the other hand, this study adopted O’Sullivan and Cui (2009) particle-particle friction coefficient and particle-wall coefficient values of 0.096 and 0.228, respectively, as the same type of spherical particles from the same manufacturer (Thomson, 2022) had been used in this study. Particle-wall friction coefficient of 0.228 was used between particles and the funnel, the set of meshes, the tube and the specimen’s mould.

The PFC manual provides a method to generate homogeneous specimens resulting in a non-zero material pressure and not including gravitation force (Potyondy, 2023). This procedure included a packing phase, which creates a grain assembly and a material-finalisation phase, which includes assigning the same material and contact properties defined in Table 2. The PFC-3D twin model is generated based on a random spatial allocation of particles while securing a target porosity. Then, the friction coefficient between particles is set to zero in order for the particles to rearrange until the mean stress nears zero or until it achieves a static-equilibrium condition. During this step, most of the overlaps between spherical particles are eliminated to prepare for the next sample preparation step. The particle size is modified until the mean stress is within the pressure tolerance of the defined sample confining pressure of 55 kPa and until the static-equilibrium condition has been achieved. In the final sample preparation phase, material properties and contact properties are defined, setting the specimen ready for shearing (Figure 7). In this sample preparation procedure, the specimen pressure is set equal to the confining pressure of 55 kPa, which is the same as the one used during the experimental triaxial test. The initial porosity is set to any random user-specified value and the target porosity is set equal to the global porosity value obtained from the experimental sample preparation and the pressure tolerance ratio is set equal to 0.1. Previous studies have also widely used radius expansion method as the sample preparation technique to simulate triaxial tests in DEM (Jenck et al., 2009; Belheine et al., 2009; Dong et al., 2020).

Since the purpose of this study is to investigate the effect of the specimen’s heterogeneity on its stress-strain response, it is hypothesised that the void ratio (⁠$e$⁠) can serve as a parameter to assess the degree of heterogeneity of the testing specimens (Hu et al., 2022; Nie et al., 2022).

In the twin PFC-3D model, a “planar” void ratio is introduced to assess the spatial variation of the specimen’s heterogeneity, by computing the void ratio along a horizontal cross section (i.e. along the x–z plane) and by repeating the same computation along different heights of the specimen (i.e. along the y direction). The variation of the “planar” void ratio along the numerical sample height can be compared with the global void ratio obtained from the sand experiment and with the PFC-3D suggested sample preparation method (i.e. radial expansion). Furthermore, the average of the “planar” void ratio along the specimen height can be considered as a representative global or volumetric void ratio in the specimen.

Before the planar void ratio is computed, the cross-sectional image on the x–z plane at different slice height, $shi$⁠, along the y axis has to be reconstructed. The distance, $di,j$⁠, between the centroid coordinates of each spherical particles above the base of the sample, $hi,j(x,y,z)$ and the slice height is calculated to check if a particular particle is within a radial distance (3 mm or 6 mm in the two numerical experiments) from $di,j$ at one particular slice height. The calculation is repeated for different slice height along the height of the specimen. The cross-section along the x–z plane of each individual spherical particles is the same for the same planar distance between the centroid of the particles to the slice regardless of its location respective to the spherical particle (i.e. the slice could fall above or below the spherical particle centroid) since the spherical particles in the twin PFC-3D model are perfectly spherical. Figure 8 illustrates the calculation process of the planar void ratio of the specimen at various slice heights. Subsequently, the coordinates of the centroid of each spherical particle were extracted and processed for calculations in MATLAB (Noble, 2012). If a particular particle is within a radial distance from $di,j$ (i.e. $di,j$< 3 mm or 6 mm), the radius of the cross section of the spherical particle is computed using the Pythagorean theorem.

In this case, $R$ is a constant value of 1.5 and 3 mm. The subscript $i$ refers to variables calculated at every slice height and $j$ refers to variables calculated for every spherical particle within the slice plane. Figure 9 depicts a vertical cross-section of a spherical particle with a horizontal slice. The sum of these projected areas is denoted as $Asi$ and the cross-sectional circular area of the triaxial mold of known diameter of 72 mm is denoted as $At$⁠. The planar void ratio, $epi$⁠, was then computed using equations (2). Binary images were also generated to visualise a slice with intersected spherical particles at various depths along the specimen.

All of the PFC simulations are run in a Windows server built on a computer equipped with an Intel Xeon E5-2695 v3 2.3 GHz processor with 14 cores (28 threads) and 128 GB of random-access memory. PFC3D version 5.0 license allows users to run two simulation instances simultaneously. The configuration of the Windows server is recommended by ITASCA personnel for running complex PFC models. It is known that computational time is an important criterion to justify if the numerical models are efficient in replicating material behaviours.

For each experimental sample, the void ratio is calculated according to equation (1). Three control variables were considered to assess the effect of sample heterogeneity on the specimen’s triaxial compression response: diameter of the spherical particles, number of diffusing sieves (placed at the top of the specimen) and particles’ drop height. The experimental design considered performing five tests under the same combination of control variable conditions. Descriptive statistics of the experimental results are summarised by the use of box plots, showing the variation of void ratio for each sample preparation condition.

Figures 10 and 11 present these summaries of the experimental descriptive statistics for the 3-mm and 6-mm diameter samples, respectively. Where the median (red line), the first and third quartiles (blue box), as well as the maximum and minimum void ratio values (black skewers) for each sample preparation condition are presented. The experimental triaxial global void ratio is also marked as bright green triangle for each experimental set up. The overall trends observed from these figures are that the sample void ratio tends to decrease as the drop height increases (regardless of the particle’s diameter), meaning that the sample becomes denser due to higher potential energy; also the void ratio of the samples prepared at the same drop height increased as more diffusing sieves are used to reconstitute the specimens, which indicates that particle’s potential energy is reduced due to contact with the sieves, leading to a lessened density of the specimen (regardless of the particle’s diameter).

One specimen for each experimental design condition (Figure 3) was subjected to a dry vacuum-consolidated triaxial compression test. Plots of deviatoric stress versus axial strain corresponding to these conditions are shown in Figures 12, 13, 14 and 15, corresponding to the 3 and 6-mm particle diameters. These figures show that the specimens behave like dense granular specimens, since all of the stress-strain curves present a clear hardening and softening behaviour, from which a peak stress can be identified.

Figures 16 presents the ultimate experimental stress of each triaxial experiment versus the corresponding specimen’s void ratio. The ultimate stresses vary from 72 to 82 kPa for 3-mm diameter sphere samples and from 67 to 85 kPa for 6-mm diameter sphere samples. Despite the fact that the void ratios of samples made of 6-mm diameter spheres are higher than the void ratios of samples made of 3-mm diameter spheres, the variation of the ultimate stress is almost similar. The effect of drop height on the ultimate stress of both specimens prepared using 3 and 6-mm diameter spherical particles, is consistent with the expectation that denser specimens are formed with higher drop height. Furthermore, the increasing number of sieves also lowered the ultimate stress of the specimen as the energy loss during contact with the sieves led to looser specimens. Nonetheless, the 6-mm samples prepared from 0.305 m drop height showed inconsistent behaviour relative to the rest of the experiments as the effect of number of sieves seemed non-existent. From the results in Figure 16, increasing the drop heights have greater influence in preparing denser specimens than decreasing the number of sieves in the sample preparation for both particle sizes. Specimen formed with 6-mm diameter particles exhibited greater increase in ultimate stress compared to the 3-mm diameter particles specimen as the sample got denser with higher drop height.

Plots of the void ratio variation of PFC twin models using the planar void ratio analysis for the numerical sample comprised of 6-mm spherical particles are presented below in Figures 17–19. The planar void ratio was computed along the height of the specimen at an interval of 1 mm. The blue dotted curve and the solid blue line represent the planar void ratio variation and its corresponding overall mean of the specimens pre-pared by the air pluviation method, respectively, while the pink curve and the solid pink line represent the planar void ratio variation and its corresponding overall mean of the specimens prepared by the radius expansion method, respectively, the most commonly used method in DEM modelling to achieve desired initial state of condition of the specimen. The mean of the “planar” void ratio can be considered as a representative value for the global void ratio of the numerical triaxial samples generated using both air pluviation and radius expansion. The green dotted curve represents the global void ratio of the specimen as obtained from the experimental observations. Consequently, these mean values from the numerical samples can then be compared to the experimental global void ratio.

As shown in Figures17–19, specimens prepared by radius expansion overestimated the planar void ratio air pluviation underestimated the planar void ratio, when compared with the experimental void ratio. All of the parameters used in the numerical experiments were set to be equal to the physical spherical particles and the geometry of the both physical and numerical experiments were set to be identical. The difference in planar void ratio is due to the inherent difference in the sample preparation technique as the radius expansion method aims to reconstruct/recreate sample given certain conditions while the air pluviation technique replicate the natural deposition process. The distribution of the void ratio greatly varies near the sample’s base due to a tight packing of the balls at the base. Moving up the height of the sample, the variation of the void ratio de-creases. At the top, the void ratio dramatically increases since there are fewer spherical particles to fill the voids (Fu and Dafalias, 2011; Hu et al., 2022). This means that the suggested sample preparation method by PFC3D is expected to produce specimens of lower strength and higher compressibility than the actual experimental behaviour, while the dry pluviation twin PFC3D method is expected produce specimens of higher strength and lower compressibility than the actual experimental behaviour.

Empirical cumulative density function (ECDF) of the planar void ratio along the specimen height for both air pluviation and radius expansion are shown in Figures 20–25 for each drop height to quantify the degree of distribution of the planar void. In all of the numerical experiments, sample prepared by air pluviation exhibited higher concentration around the mean value of the planar void ratio as shown by the maximum gradient of the ECDF curves. However, due to the significant variation in the planar void ratio at the bottom and the top of the specimen, there is negligible difference in the overall sample standard deviation between samples prepared by the two preparation methods. Furthermore, the mean of the planar void ratio for each sample preparation method per numerical experiments are compared to the experimental void ratio in Figure 26. To quantify the difference between the two sample preparation techniques and the experimental void ratio, the absolute percentage error is also computed by taking the difference of the void ratio values and dividing it by the experimental void ratio. The results are also presented in Figure 26. Together with the distribution shown from the ECDF plots, it can be concluded that sample preparation by air pluviation yields sample that has a more uniform planar void ratio and closer to the experimental void ratio. Consequently, the DEM model using air pluviation method is expected to provide more accurate sample responses when simulating compression test.

Computationally, samples prepared in PFC by the radius expansion method take a couple hundred thousand cycles, approximately 1 h, to complete while samples prepared by the air pluviation method take 30 m cycles, approximately three days, to complete. The difference is due to the fact that the air pluviation method simulates the actual physical stage of each sphere over time thus consuming more calculation cycles. The calculation time of the radius expansion method can be increased or decreased by changing the pressure tolerance.

This study highlights the sample mechanical behaviour dependence on the sample preparation method. Three control variables are defined to assess such dependence: diameter of spherical particles, number of diffusing sieves and particle drop height used to reconstitute specimens.

From a micro-scale perspective, the analysis methodology presented in this research article presents a systematic assessment of the microstructure characteristics of a particle-based specimen prepared by various methods (e.g. air pluviation and radius expansion). Analysing the microstructure of the resulting specimen can further our understanding of the controlling macro behaviour of a particle-based materials

This paper shows a significant difference when using the radius expansion method versus a sedimentation-based process like air pluviation, which raises serious concerns about likely risk impacts associated to using what the literature shows as a “default” method for DEM model preparation, because of the known effects of voids in terms of materials strength, compressibility and permeability. Consequently, results presented in this paper contribute to reduce risk impacts when using DEM models to simulate macro-scale structures that pose a significant threat to human life, the environment and the economy (e.g. foundations of civil infrastructure in static and dynamic conditions, landslides, earth dams, tunnels, etc.), since it anticipated a potential underestimation of materials strength, overestimation of materials compressibility and overestimation of materials permeability (which in geomaterials also affects negatively strength and compressibility). Moreover, the presented arguments show a systematic deviation of voids estimates when comparing radius expansion versus air pluviation, reaching up to twice the error. Furthermore, this study highlights that the reliability and quality of data generated through DEM are strongly influenced by the chosen modelling approach, underscoring the importance of carefully selecting methodologies that yield physically representative data sets, that can be used later to train surrogate models.

One of the key advantages of DEM lies in its capacity to capture large deformations and localised failure mechanisms that continuum-based methods often cannot. This capability positions DEM as a powerful and efficient tool for generating comprehensive databases and realistic representations of failure mechanisms. When initial conditions are defined in closer alignment with physical reality (e.g. shear-bound conditions), data-driven methods can more effectively learn from these simulations, yielding improved predictive models. In turn, this facilitates the exploration and accurate capture of both global and local effects while enabling the generation of extensive training data sets without the need for costly and time-intensive physical experimentation.

The corresponding author would like to acknowledge the funding provided by the Zachry Career Development Professorship II at Texas A&M University to complete this manuscript