3D printed porous particle and its geotechnical properties

The study of effects of particle breakage on the mechanical properties of soil composed of porous particles is challenging due to the heterogeneity of the shape and inner void structure of individual particles, even for an identical soil sample, which imparts a compound effect on the mechanical properties. Advancements in three-dimensional (3D) printing technique have enabled the replication of objects with the same shape but different inner structures. This study investigated the feasibility of replicating porous and non-porous particles with the same particle shape characteristics, such as form, waviness, and texture, using 3D printing technique. The particle shape characteristics were evaluated using image analysis. Single particle crushing and triaxial compression tests were conducted to characterize the mechanical properties of the 3D printed and porous volcanic soil particles. It is observed that the mechanical response in the single particle crushing test varies for volcanic soil, which may be attributed to the heterogeneity in the shape and porosity of the particles. However, for each type of 3D printed particle, the response has a high repeatability and varies based on particle porosity. Furthermore, the effects of porosity on the shear response are demonstrated through triaxial tests on 3D printed particles of different porosities. It is noted that although a quantitative comparison is not possible, a qualitative similarity is observed in the response of the 3D printed porous particles with natural porous volcanic soil. Thus, insights into the mechanical response of porous particles can be gained using 3D printed particles.


Introduction
Mechanical properties of granular materials, such as soil, are characterized by various factors, including packing density, degree of saturation, stress condition, gradation, shape, and strength of constituent particles.These factors influence each other to characterize the mechanical properties.Therefore, the evaluation of the individual effects of particle shape and its strength is challenging, unlike that of other factors that are controllable.This is because natural soils are composed of irregularly shaped particles with a variety of strengths resulting from mineral composition or inner structures such as pores.Generally, it is impossible to test multiple specimens composed of particles of the same shape and strength, even if these specimens can be tested at the same gradation, packing density, degree of saturation, and stress condition.This is because tests to characterize the mechanical properties of soil, such as triaxial compression tests, lead to an irreversible destruction of the particles.
Numerical simulations using the discrete element method (DEM) have significantly contributed to solving this challenge.In DEM, irregularly shaped particles can be modeled using aggregates composed of multiple smaller spheres [8,50] or based on particle surface geometry [23,49], where these strengths can be considered as an input parameter.Although the shape is simplified to some degree, depending on the packing method of the smaller particles or the resolution, to represent particle surface [9,22,44], the typical shear response of soil has been successfully simulated.Moreover, aggregates modeled by bonding smaller spheres enable the simulation of particle crushing behavior [8,11,13,41,46,50].As investigated by Xu et al. [46], modeling of breakable and unbreakable particles with the same particle shape and different bond parameters helps us understand the net effect of particle breakage.Wang et al. [41] proposed a calculation method to quantitively evaluate the degree of particle breakage using a spherical aggregate to directly compare a corresponding parameter in experiment called breakage ratio.Although the crushable particles can be modeled and their behavior can be simulated by DEM, a systematic experimental study is needed to verify results of the simulation and to develop more accurate force prediction models.
In recent years, three-dimensional (3D) printing techniques have rapidly advanced and achieved a high printing resolution with a variety of printing materials.Accordingly, the technique has been used to evaluate the properties of soil or rock using 3D printed particles or pore structures.Several studies have replicated pore networks of loam for hosting fungus [36] and for hydraulic conductivity [7].Hanaor et al. [14] investigated the feasibility of 3D printed irregularly shaped particles modeled systematically using few mathematical functions to verify whether the typical shear response of particles is similar to that of granular materials.Researchers have revealed various properties of 3D printed soil or rock analogs related to shear behavior [1-3, 15, 31, 38], hydraulic conductivity [43], wettability of printing material [20], and transparency of 3D printed particles [29].In addition, 3D printing technique is expected to be used for calibration of DEM simulation [11,13,26,30,37] by modeling the same particles through experiments and simulations.A significant advantage of the 3D printing technique is that particles of the same shape but with different inner structures can be replicated.Ge et al. [10][11][12][13] investigated the feasibility of 3D printed porous particles for particle crushing behavior using a single particle compression test.Xu et al. [47] revealed the effects of inner fissure on creep behavior during 1D compression using 3D printed irregularly shaped particles.However, in the literature, the individual effects of inner pores of particles on the shear response or a systematic methodology to create the 3D printed particles with the same shape but different inner structures were not clarified.
Hokkaido, Japan, is a well-known land widely covered by volcanic soils, and such soils have triggered serious geodisasters induced by earthquakes [19,24,25,42].Coarse volcanic soil particles are characterized as highly porous with low particle strength.The literature [17,18,33,35] suggests that the porosity of volcanic sand is a more concerning parameter that can significantly affect the mechanical response than other parameters.Its evaluation further becomes complex as heterogeneity in the material results in significant variation in the particle strength.In general, several soil tests are required to characterize the mechanical properties for each type of volcanic soil.Therefore, an experimental method to evaluate the porosity or particle strength as a comparable or controllable parameter in a systematic manner is required.
The study of the effects of particle breakage on the mechanical properties of soil composed of porous particles, such as volcanic soils, remains challenging.Conversely, 3D printing can be a promising method to replicate synthetic porous particles and enable tests to verify the effects of particle breakage in a more systematic manner.However, there are very limited studies in this field [10][11][12][13]47].Moreover, no methodology for synthetic porous particle modeling has been established, and mechanical properties such as shear response have not been investigated.This study aims to assess the feasibility of 3D printing techniques to create particles of the same shape but different inner structures (i.e., voids) and investigate whether these particles exhibit the typical mechanical behavior of porous particles, such as crushable volcanic soils.Therefore, we created various types of 3D printed particles modeled using DEM simulations and conducted single particle crushing and triaxial compression tests in addition to a natural porous volcanic soil (VS).The three types of 3D printed particles are as follows: sphere particles (3DPs), particles with inner voids (3DPv), and particles without inner voids (3DPf).3DPv and 3DPf were modeled and printed ensuring they had the same shape but different inner structures.Subsequently, we described the steps to model 3D printed particles based on the DEM simulation and their fundamental physical properties and those of VS.Moreover, we examined the image analysis results of the particles to verify the formation of the intended particles.Finally, we examined the results of the mechanical tests and particle behavior.
Figure 1 presents the procedure to model the three types of 3D printed particles: 3DPs, 3DPv, and 3DPf.The models were generated using the commercial DEM software Particle Flow Code (PFC) ver.5.0.First, a spherical particle with a diameter D 1 was generated at Step 1, corresponding to 3DPs.The second step was to pack smaller unit spheres with different diameters (D 2 and D 3 ) at the target porosity and volume fractions in the spherical particle of diameter D 1 and to allow them to overlap each other (Step 2).At this step, the sphere particle with D 1 is used only as the boundary where the unit spheres can be filled, and therefore it does not actually exist for 3DPv and 3DPf after the step.The unit spheres of diameter D 3 were subsequently removed, and the aggregate composed of unit spheres of diameter D 2 remained (Step 3).The aggregate generated at Step 3 can be a model of 3DPv that has inner voids, which is created by removing the unit spheres of diameter D 3 in addition to the original inter-particle voids that existed at Step 2. Notably, because the unit spheres forming 3DPv were allowed to overlap, the aggregate did not loosen up and maintained its shape, even after it was 3D printed.Subsequently at Step 4, 3DPf was modeled by generating a single unit sphere with a diameter D 4 at the center, which was slightly smaller than D 1 , on the 3DPv model.This process resulted in inner voids in 3DPv filled by the unit sphere of diameter D 4 , while some unit spheres with D 2 protrude from the one with D 4 as shown in the figure.Consequently, we assume that 3DPv and 3DPf have identical particle shape characteristics with different inner structures.Moreover, 3DPv and 3DPf seem to have the same particle size as 3DPs, because these were generated by filling 3DPs with the smaller unit spheres.The similarity in the particle shape characteristics was verified using the image analysis discussed in the subsequent sections.
Table 1 lists the parameters used to model the 3D printed particles in PFC.The parameters D 1 -D 4 are the diameters used to generate the unit spheres, as indicated in Fig. 1.D 1 is 5 mm, which is the diameter of 3DPs and is also assumed to be equal to the circumcircle diameter of 3DPf and 3DPv.The sphere with diameter D 1 was randomly filled with unit spheres of diameters D 2 (0.6 mm) and D 3 (0.4 mm), with individual volume fractions of 0.65 and 0.35, respectively, to achieve a porosity of 0.36.These values were determined by trial-and-error modeling and 3D printing such that 3DPv could maintain its aggregate without the unit spheres being unglued after 3D printing before performing the subsequent tests.D 4 (4.4 mm) was determined as the difference of D 1 and D 2 .This indicates that the unit spheres of diameter D 2 located at the outer rim of 3DPf exits the unit sphere of diameter D 4 at a maximum distance of D 2 /2.This is attributed to the fact that a significantly large D 4 appears to start making contacts between the unit spheres of diameters D 2 and D 4 during the triaxial compression test; these contacts are undesired because 3DPf must have a shape and skeletal structure similar to that of 3DPv.On the contrary, if D 4 is considerably small, the outer unit spheres can be easily unglued, resulting in particle shape characteristics different from those of 3DPv.Table 2 presents the parameters of the constituent spheres for the 3D printed particles: the diameters (D 1 -D 4 ) and numbers of unit spheres with diameters D 1 -D 4 (n 1 -n 4 ).As a result of particle generation, n 2 and n 3 Fig. 1 Manufacturing process of 3D printed particle models Acta Geotechnica (2023) 18:5735-5753 5737 were 268 and 497, respectively, whereas both n 1 and n 4 were 1.
Because 3DPv has been modeled using several unit spheres, the strength of a single 3DPv is possibly dependent on the strength to bond adjacent unit spheres.Ge et al. [10] introduced aggregates modeled on regularly arranged unit spheres and the bonds bridging them.Consequently, the bond strength can be controlled using the bond length (BL) (Fig. 2a) and material properties that can be used to individually print the bonds.Moreover, the BL controls the porosity of the aggregate.Particularly, the greater the BL, the higher the porosity of the aggregate.In this study, because the bonds between the adjacent unit spheres are created by overlapping spheres, the bond strength is possibly related to the overlap length (OL) (Fig. 2b). Figure 2c shows the frequency distribution of the OL for 3DPv.The number of contacts where adjacent unit spheres overlap is 522, whereas the number of unit spheres (i.e., n 2 ) is 268 (Table 2).Figure 2c indicates that several overlaps have an OL of 0.05 mm or smaller, which corresponds to lengths lower than 10% of D 2 .The percentage of OL smaller than 0.05 mm is 70.1%, and the average OL is 0.038 mm.
Each 3D printed particle model was transformed to a printable standard triangulated language (STL) model using a PFC function.Su et al. [38] reported that the extent of particle morphology simplification (i.e., the number of triangular faces to form particles) could affect the shear response of 3D printed synthetic specimens.In this study, the parameter related to the resolution of the STL model was uniformly set to 0.5 (Table 1).Therefore, every sphere including the unit spheres was modeled with the same number of triangular faces.The properties of the STL models are listed in Table 3.A single sphere, 3DPs, was modeled with 320 triangular faces; for 3DPv composed of   Including the vertices and triangular faces inside the unit sphere of diameter D 2 268-unit spheres (Table 2), the number of triangular faces was 85,760.Figure 3 shows the models represented by the STL format and printed particles.The 3D printer used was Agilista-3100 (Keyence Corporation).The manufacturer specifications report that 3D objects can be printed with a resolution of 40 lm 9 63.5 lm in the horizontal plane with a thickness of 15 lm.The 3D printed objects were created using two different UV-cured resin-based materials identified as ''AR-M2'' and ''AR-S1'' to form the solid matrix of the object and temporary support material in the object spaces, respectively.The properties of printing material (i.e., AR-M2) are listed in Table 4, and more details regarding printed objects have been reported by Matsumura et al. [31].

Volcanic soil
Natural porous VS (particle size 2-4.75 mm) was used in this study.The soil was classified as poorly graded sand (SP) as per the Unified Soil Classification system.The soil was obtained by scalping specific size fractions from a fall pyroclastic flow deposit that erupted from Mt. Tarumae in Hokkaido and is referred to as Ta-d [21,24,25,28].Ta-d has been reported as the volcanic soil that caused numerous slope failures due to the 2018 Hokkaido Eastern Iburi earthquake [24,25].and 4d) and low particle strength is illustrated later.
Table 5 lists the physical properties of the original Ta-d sampled at depths of 0.8 to 1.8 m in this study.As shown, the physical properties of the original Ta-d exhibit significant variations, even though they are sampled with a layer thickness of 1.0 m.This may be attributed to the difference of degree of weathering in the Ta-d layer; that is, the properties such as particle size and mineral composition have changes depending on the depth in the layer through the weathering process [21,25,28].VS was collected from the original Ta-d with the same degree of weathering, which can be visually confirmed, and prepared by ovendrying and sieving it to 2-4.75 mm (Fig. 4c).The triaxial compression properties of the original Ta-d reported by Kawamura [24] are also presented in Table 5.The cohesion, c d , is 33 kPa and the internal friction angle, / d , is 29.0°.The increment of fines content, DF c , is 3.9-8.6%obtained based on the sieve analysis after the triaxial compression tests (not reported in Kawamura [24]).The friability of the original Ta-d can be considered moderate  Poisson's ratio 0.37 b a From Matsumura et al. [31] b From 3D printer manufacturer's literature based on the Japanese Industrial Standard (JIS), JIS K7161, Plastic-Deformation of tensile properties and in conformity to ISO527.E = 1820 MPa based on the unconfined compressoin test results from Matsumura et al. [31] compared to those of other volcanic soils sampled in Hokkaido [33].

Particle properties of the tested materials
Table 6 lists the particle properties of the four types of materials.The particle density of the 3D printed particles is 1105 kg/m 3 [31], while that of VS is 2676 kg/m 3 .VS tested for the particle density was first sufficiently ground on a mill.Subsequently, the ground VS particles were tested based on JIS A1202, where the soil was boiled in a pycnometer for more than 2 h to remove the air bubbles inside the particles.Therefore, for VS, the result in Table 6 indicates the true density of particle not including volume of inner voids.The particle size of the 3D printed particles refers to its diameter, D 1 as described in Sect.2.1, whereas  Based on the sieve opening size that of VS is based on the sieve opening size.Adamidis et al. [1] reported that the particles should be sized up (more than 2 mm) and printed accordingly to capture any irregular particle shape characteristics and mechanical properties of Hostun sands smaller than 1 mm.Based on their suggestion, 3DPs, 3DPf, and 3DPv are sufficiently large in terms of their diameters.However, because 3DPf and 3DPv comprise considerably smaller unit spheres (Fig. 1), a mass analysis is performed to verify whether the resulting printed particles have different inner structures.The mass of particles was calculated using the mean value of 30 particles for each 3D printed particle model and 300 particles for VS.Thus, the mass-per-particle values of 3DPs, 3DPf, and 3DPv were 0.068, 0.052, and 0.024 g, respectively.Because the mass of a sphere with a diameter of 5 mm and density of 1105 kg/m 3 theoretically must be approximately 0.072 g, 3DPs was printed with a slightly lower mass.This may be attributed to the resolution of the STL model to form the sphere (Fig. 3a and Table 3).The mass of 3DPf is smaller than that of 3DPs owing to the difference in particle asperity.The mass of 3DPv is smaller than that of 3DPf because inner voids exist in 3DPv.Thus, clearly, 3DPf and 3DPv have different inner structures based on the results of mass analysis.The porosity of particle (e) is defined as: where V s and V p are the volumes of solid part and pore in the particle, respectively.The porosity of 3DPs and 3DPf is assumed to be 0% based on the modeling method.If we assume that the difference between 3DPv and 3DPf lies only with respect to the inner structure, whereas their particle shapes are identical, V p for 3DPv (V pv ) can be assumed as: where V sf and V sv are V s for 3DPf and 3DPv, respectively.From Eq. 2, (V s ?V p ) in Eq. 1 for 3DPv is equal to V sf .Therefore, the porosity of 3DPv (e V ) can be represented as Eq. 3 and calculated by the mass ratio of 3DPv to 3DPf using the particle density and mass in Table 6 as Eq. 4.
where m sv and m sf are the mass of 3DPv and 3DPf, respectively.Eventually, the porosity of 3DPv is calculated as 53.8% from Eq. 4 (Table 6).
The porosity of VS (59.9%) was measured by the gas replacement method using AccuPyc II 1345 and beads replacement method using GeoPyc 1365, both of which are manufactured by Micromeritics Instrument Corporation.Those methods enable the measurement of V p and V s for VS, and its porosity is calculated based on Eq. 1 using the measured volumes.The tests were repeated five times using 5-8 g of VS dried at 110 °C and vacuum-degassed for 24 h and the porosity listed in Table 6 was calculated as an average of the results.Notably, although the porosities of 3DPv and VS are relatively close, their pore size distributions are considered significantly different.The SEM image (Fig. 4d) indicates that most of the voids for VS are smaller than 0.05 mm.According to the results of the mercury intrusion porosimetry (not reported in the paper), the voids of a VS particle occupied more than 50% with pore diameters of 0.001-0.010mm.In contrast, because 3DPv was modeled by removing the unit sphere of diameter D 3 inside the aggregate (Fig. 1), its pores were possibly larger than the diameter of the removed unit spheres, that is, 0.4 mm (Table 1).
3 Testing methods

Image analysis to investigate particle shape characteristics
The effect of particle shape on the mechanical properties of granular materials is well-established in literature [4-6, 39, 45].The image-based particle shape characterization method proposed by Vangla et al. [39] was used to ensure the similar shape of the 3D printed particles used in this study.To this end, digital microscopic images of three particles were captured from three directions using a digital microscope system, including VHX-6000 and VHX-S650 manufactured by Keyence Corporation.Each image has 1600 9 1200 pixels.The images were converted to binary images, and operations were performed on it using a MATLAB code to quantify the parameters.In this study, the different parameters, such as circumcircle, sphericity, roundness, and normalized roughness were used to quantify shape of the particles at different scales (macro, meso, and micro) [39] (Fig. 5).Sphericity is defined as the widthto-length ratio sphericity used in the visual chart proposed by Krumbein and Sloss [27].The roundness is computed Acta Geotechnica (2023) 18:5735-5753 5741 using the formula proposed by Wadell [40] and is defined as the ratio of the average radius of all convex regions of the particle shape that are less than or equal to the radius of the maximum inscribed circle.The normalized roughness parameter proposed by Vangla et al. [39] is the root mean square roughness of the particle normalized by the length of the particle.

Single particle crushing test
A series of single particle crushing tests was performed using a typical uniaxial compression test apparatus (Fig. 6) to investigate the particle crushing strength of the porous material.An air-dried particle was placed in a stable direction on the bottom platen, and the top platen was allowed to contact the particle.A vertical load was applied by moving the bottom plate at a rate of 1 mm/min until particle breakage was visually confirmed or a vertical displacement of approximately 2 mm was reached.A digital camera was used to capture some particle images from the side during loading for visual inspection and presentation.A total of twenty particles for each type of material was tested.Every test was conducted using the newly generated 3D printed particles and fresh VS.

Triaxial compression test
A triaxial compression test was performed to characterize the mechanical properties of the four types of materials.Samples with 50 mm in diameter and 100 mm in height were prepared to conduct the tests.Identical connections and skeletal stresses exerted through the contacts were ensured by preparing each specimen of 3D printed particles by packing the same number of particles into each specimen.Table 7 presents the testing conditions for all the samples.For 3DPs and VS, the specimens were packed at minimum dry densities by depositing particles using a funnel with an opening diameter of 20 mm.The number of particles packed in the specimens was assumed to be 1854 for 3DPs and approximately 4000 for VS based on the measured particle mass (Table 6).For 3DPf and 3DPv, the specimens were prepared by compaction to contain the same number of particles as that of 3DPs.Compaction was necessary because the particle asperity of 3DPf and 3DPv prevented densification to some extent, unlike 3DPs.However, compaction was performed with a low compaction energy so as not to break or abrade the particles.All particle types were packed after being air-dried.After packing the particles, the specimens were fully saturated to achieve a B value higher than 0.95.Subsequently, the specimens were isotropically consolidated at target stresses of 50, 100, or 200 kPa.The consolidation process continued for 24 h.After consolidation, an axial loading with a rate of 0.5%/min was applied under a drained condition and constant isotropic confining stress until the axial strain reached 15%.During axial loading, the axial displacement and axial load were monitored using a displacement transducer and load cell installed out of the cell.Conversely, the volume change was monitored using a burette connected to the top and bottom of the specimen.Fresh particles were used for every porous 3DPv and VS test because some of the tested particles showed particle breakage, whereas the particles of 3DPs and 3DPf, which is non-porous, were partially reused.However, tests were performed from lower to higher confining stress conditions to minimize the reuse effect.Adamidis et al. [1] reported that the effect of reuse of 3D printed particles was insignificant; however, owing to the limitations of the study, more investigations with different types of particles and stress boundaries are necessary.

Particle shape characteristics
Figure 7 presents the shape descriptors obtained for the different particles using image processing.Figure 7a shows  the circumcircle diameter for all samples.Evidently, the circumcircle diameters of 3D printed particles are comparable to those of VS.Furthermore, we observed that the variation, depicted by the sample standard deviation, was insignificant for the 3D printed particles compared to that for VS.The sphericity for the 3D printed particles is shown in Fig. 7b; 3D printed particles have similar sphericity of approximately 1.0.However, the average sphericity of VS is approximately 0.75, and a higher variation is observed in comparison to that in the 3D printed particles.Figure 7c presents the normalized roughness of the samples; the normalized roughness for the 3D printed particles is within a similar range, and it is less than that of the natural particles (i.e., VS). Figure 7d shows the roundness of the samples; evidently, all 3D printed particles have roundness within a similar range, whereas VS has a comparatively higher roundness.Notably, the variation in the roundness value of 3DPs is the highest among all the 3D printed particles, and it is higher in comparison to the variation in all the other descriptors for all the 3D printed particles.This may be attributed to the complications involved in 3D printing of the spherical particles.However, considering all the variations and complications involved, the descriptor values validate the assumption of printing particles with different internal structures and similar shape.Thus, the 3D printed particles have similar shape characteristics.However, the particle shape characteristics of VS, which is a natural soil and irregularly shaped, are different from those of any 3D printed particles and have higher variation.

Crushing behavior of single particles
A series of single particle crushing tests were conducted on different particles selected from VS to understand the effect of porosity on the mechanical response of natural porous soil.Figure 8 presents their axial force vs. axial displacement responses.It is evident from Fig. 8a that the force-displacement responses and magnitude significantly vary for VS owing to the variation in porosity and particle shape.The maximum force value varies widely in the range of 5 N to 60 N, and the secant modulus at 0.5 mm axial displacement ranges from 3 to 55 N/mm.However, these responses can be broadly classified into three categories: (i) an initial increase in the axial force followed by strain softening and hardening (Fig. 8b), (ii) a linear increase in the axial force with an increase in the axial displacement (Fig. 8c), and (iii) an initial linear increase followed by strain hardening (Fig. 8d). Figure 8b presents the trendline fitted (polynomial curve of order 3) to the force and axial displacement responses of category (i) to demonstrate the general trend.The porosity could be responsible for these different categories of responses of VS; however, it is not explicit as the tested particles have different shapes.Thus, a series of experiments was conducted on three types of 3D printed particles with different porosities and the same shape to systematically understand the effect of porosity on the mechanical response of porous soil.Figure 9 presents the axial force vs. axial displacement responses of the 3D printed particles.Figure 10 presents the secant modulus at 0.125 mm (E 0.125 ), 0.25 mm (E 0.25 ), 0.5 mm (E 0.5 ), and 1.0 mm (E 1.0 ) axial displacement for 3D printed and VS particles, respectively.The average secant modulus with deviation is provided for the 3D printed particles; however, the range bound of secant modulus is provided for VS to demonstrate the variation in the values due to the heterogeneity of the natural particles.Figure 9 shows the repeatability in the responses for each 3D printed particle, which is unlike that of VS, representing the potential of obtaining 3D printing particles with identical shapes and porosity.The response varies with changes in the porosity of the 3D printed particles.displacements.The axial force initially increases linearly with an increase in axial displacement until the lower yield is subsequently followed by yield strength.Subsequently, a plateau of constant force is observed due to yielding until the strain hardening/increase of axial force commences, which can be attributed to the densification of the 3D Fig. 8 Axial force and axial displacement responses of the single particle crushing test for VS: a all tested particles; the general response of particles showing b an evident peak followed by softening and hardening, c linear increase, and d linear increase followed by strain hardening printed particle matrix [16] during the yielding stage and increase in contacts between the particle and load platen.A similar response is observed in category (i) of VS particles.It is observed in 3DPf initially that the increase in the axial force rate is low compared to that in 3DPs, which is also manifested by the secant modulus (E 0.125 and E 0.25 ) presented in Fig. 10.However, as the axial displacement increases, the difference in the secant modulus values decreases, as observed for E 0.5 and E 1.0 .This is attributed to the initial compression of outer unit spheres of 3DPf or their breakage during the axial compression.Further, it is noted that the strain hardening beyond yielding is comparatively higher in the case of 3DPf, manifested from the rate of change of slope [39] beyond 1.5 mm displacement.For instance, the average rates of change of slope are 301.2 and 287.3 N/mm for 3DPf and 3DPs, respectively.The higher strain hardening in 3DPf is due to the densification of the structure after the initial compression of the asperities on the surface.Figure 9c shows that for 3DPv, the axial force increases linearly up to 1.5 mm displacement, unlike for 3DPs and 3DPf, where a higher rate of increase in axial force is observed.Furthermore, the magnitude of the force is approximately 30 times less than in the case of 3DPs and 3DPf.The secant modulus determined at 0.125, 0.25, 0.5, and 1 mm (Fig. 10) are almost similar for 3DPv, representing an initial linear response, and the values are significantly lower in comparison to those of 3DPs and 3DPf due to high porosity.However, around 1.5 mm, the  In comparison to VS, both 3DPs and 3DPf exhibit a higher yield strength, whereas 3DPv exhibits a low yield strength in the range of 5 N to 15 N, which is comparable to that of some VS particles.Unlike in VS, no sudden drop in the axial force is observed in the case of 3D printed particles.This is attributable to the ductile nature of the resin resulting in gradual failure of the particles [47].However, a similarity in the response for different 3D printed particles is observed with the categories classified for VS.For instance, 3DPs and 3DPf show a linear increase followed by strain softening and hardening like the first group of VS, as shown in Fig. 8b.Notably, 3DPv shows a linear increase followed by strain hardening as shown by some particles of VS in Fig. 8d.Furthermore, the magnitude of 3DPv is very low compared to those of 3DPs and 3DPf as also observed in VS with a similar response and low magnitude compared to other particles of VS.It is observed from Fig. 10 that the secant modulus values at different displacements vary widely for all the VS particles, unlike 3D printed particles.The secant modulus values are comparatively lower than that of 3DPs and 3DPf for all the VS particles.However, a similar trend of initial increase and decrease in secant modulus is observed for VS particles of the category (i), 3DPs, and 3DPf.Further, the values of secant modulus for the category (ii) and (iii) of VS particles are similar to 3DPv representing a linear response.In conclusion, although the 3D printed particles are ductile and difficult to break, a similarity in the response is observed for both 3D printed and VS particles with change in the porosity, and thus, insights into the effect of porosity on the mechanical response of natural porous volcanic soils can be gained using 3D printed particles.
Figure 11 shows images of each particle captured before and after loading to examine the failure mode.Notably, VS placed on a table was horizontally compressed using a caliper only to examine the failure mode, because VS cannot retain crushed pieces on the testing apparatus during the test.It is visually evident from the figure that there is no particle crushing for the 3D printed particles (Fig. 11a, b,  and c), even though the top and bottom of the particles contacting the platens are compressed to a flat state.The horizontal displacement measured at the most inflated state, which is calculated from the images, is approximately 0.07 mm for 3DPs, 0.05 mm for 3DPf, and 0.12 mm for 3DPv.3DPv horizontally expanded more under the vertical load.Therefore, the tensile fracture by horizontally splitting the particle appears to occur inside the particle for 3DPv, although it is not apparent in Fig. 11c.The behavior is considered to cause a significantly lower particle strength (Fig. 9c).Although particle breakage was visually observed for all tested particles of VS, the failure pattern differed according to individual particles-few particles showed compression failure and were crushed to several pieces (Fig. 11d), whereas others exhibited a tensile failure and were split into a few pieces (Fig. 11e).The tendency likely depends on particle shape characteristics, boundary conditions, and the inner void structure.

Triaxial compression behavior
Figure 12 shows volumetric strain behavior during an isotropic consolidation.It is observed from Fig. 12a that the increase in the contractive volumetric strain with the increase in confining stress for 3DPs is insignificant.The overall volumetric strain is less than 5% at all confining stresses for 3DPs.In comparison to 3DPs, 3DPf exhibits higher contractive volumetric strain (Fig. 12b), and also the difference in the volumetric strain slightly increases as the confining stress increases.The higher volumetric strain in 3DPf compared to 3DPs may be attributed to the asperity breakage related to the plastic behavior or rearrangement caused by slipping at the contacts of unit spheres.For porous 3DPv, the contractive volumetric strain is more significant than that of 3DPf, even at the lowest confining stress of 50 kPa.The contractive volumetric strain attains a value of approximately 30% at 200 kPa.Thus, the results indicate a significant effect of particle porosity on the contractive behavior due to isotropic consolidation.For VS, the volumetric strain is significantly low at 50 kPa compared to those of the 3D printed particles; however, it increases considerably as the confining stress increases.The results suggest that the stress-dependent tendency is common for porous particles (i.e., 3DPv and VS).An important difference between the 3D printed particles and VS is that the volumetric change of the 3D printed particles continues toward the end of consolidation, whereas that of VS converges to indicate the end of consolidation.Such time-dependent behavior of the 3D printed particles was also experimentally confirmed [2,47].A possible reason for the behavior is related to more plastic inter-particle contact response due to the deformation of particle asperities for the 3D printed particles [2] in addition to a high contraction potential of the porous particles [47].
Figure 13 shows the relationships of stress ratio calculated by r 1 /r 3 and volumetric strain with axial strain during shear loading.For 3DPs and 3DPf, the peak strengths appear to exhibit a stress-softening behavior for all confining stress conditions.Notably, the peak and residual strengths of 3DPs are apparently lower than that of 3DPf.The difference in particle asperity resulted in a difference in these strengths.The results of the volumetric strain indicate that the dilative behavior is observed after a slight contractive behavior at the beginning of loading.The effect of confining stress is insignificant for 3DPs, while it appears for 3DPf, especially between 50 and 100 kPa.Further, it is observed that the initial stress ratio (before 6% in axial strain) in 3DPs is higher than 3DPv; however, beyond that the stress ratio is higher in the 3DPv.The lower stress ratio values in 3DPv than 3DPs until 6% in axial strain may be attributed to the crushing of the 3D pores.However, with an increase in the deviatoric load, the particles may have undergone a densification process, resulting in the microstructure rearrangement in 3DPv particles and, thus, higher shear strength.Further, higher stress ratio values in 3DPv beyond 6% in axial strain may also be attributed to the asperities on their surfaces.The effect of porosity is observed in the shear response of 3D printed particles.It is observed that with the increase in porosity, the stress-softening has reduced, which is also observed in VS.A progressive failure is observed for the porous 3D printed particles (3DPv), which is similar to natural porous particles (VS).Further, the volume change is only compression in porous materials, 3DPv and VS, unlike 3DPs and 3DPf, where the volume changes from compression to dilation.Moreover, the porous 3DPv and VS exhibit a decrease in strength as the confining stress increases.Such stress hardening and stress dependency are typical behaviors for crushable soils including volcanic soils [17,18,33,48], and are also observed in DEM simulation [46].Xu et al. [46] modeled breakable and unbreakable agglomerates in hexagonal close packing with the same particle shape and different bond strengths and investigated the net effect of particle breakage on the shear behavior.Interestingly but not surprisingly, the breakable agglomerate exhibited the above-mentioned stress dependency and contractive dilatancy behavior like 3DPv and VS, although the stress level is different (i.e., it is lower in this study).As seen in Fig. 13c and d, an interesting difference is that the amount of volumetric strain for 3DPv is significantly lower than that for VS.This was possibly because a more contractive volume change had already occurred during isotropic consolidation for 3DPv compared with that for VS, as shown in Fig. 12.Thus, there was little scope for further contraction during shear loading.In contrast, for some coarse-grained volcanic soils, more particle breakage tends to occur during shear loading than during isotropic consolidation [18].Particle breakage generally leads to volume contraction.Therefore, more contraction occurs during shear loading for VS.
Based on the results of triaxial compression test with three different confining stress for each material, the strength parameters--that is, the cohesion (c d ) and internal friction angle (/ d )--are compared, as shown in Fig. 14.Evidently, a typical tendency for crushable soils appears in the result of VS; c d is high and / d is low.The particle breakage with volume contraction disturbs the increment of shear strength as the consolidation stress increases.Consequently, c d apparently increases, while / d reduces.
Comparing VS with the 3D printed particles, a relatively Fig. 12 Volumetric strain during isotropic consolidation of a 3DPs, b 3DPf, c 3DPv, and d VS higher c d was observed for 3DPf and 3DPv, while it was zero for 3DPs.A possible reason for this difference is that protrusions on the outline of particles lead to interlocking between particles and higher strength, even at lower confining stress.Moreover, in the case of 3DPv, as the confining stress became higher, more volume contraction occurred and disturbed the increment of shear strength.This may explain the lower value of c d and / d for 3DPv than that for 3DPf.
As an index to evaluate the stress dependency of crushable soils such as volcanic soils, the internal friction angle (/ 0 ) calculated for each test result assumed as cohesion = 0 is compared [33].Figure 15 shows the relationship between / 0 and the mean effective principal stress at failure, that is, at q max (¼ r 1 À r 3 ), p f .In the figure, the results of triaxial compression test including the consolidated drained (CD) and undrained (CU) conditions for few volcanic soils sampled in Hokkaido [33] and the original Ta-d [24] are presented.Evidently, depending on the type of volcanic soil, the tendency for decrease in / 0 with an increase in p f significantly differs.For example, the Nakashibetsu volcanic soil shows a higher stress dependency, in which / 0 significantly decreases with increasing p f , whereas Mori does not.It is reported that this tendency is related to the extent of particle crushing; the greater the particle crushing, the higher the stress dependency [33].However, the original Ta-d, which has moderate friability (Table 5), shows higher stress dependency.The reason for this observation is clearly explained in literature; however, the mechanical properties of volcanic soils can be considered to be highly dependent on the volcano (i.e., the eruption source).The porous 3D printed particles, 3DPv, show stress dependency especially at the lower stress level.3DPf also shows stress dependency, although the tendency is slightly lower than that of 3DPv.In contrast, no stress dependency is observed for 3DPs in the tested stress level.

Discussions on challenges and prospects for the geotechnical usage of 3D printed porous particles
The aforementioned results demonstrate that the 3D printed porous particles 3DPv have qualitatively similar mechanical properties compared to the natural porous VS in terms of lower particle strength (Figs. 8 and 9), highly contractive behavior (Fig. 12), and stress-hardening behavior (Fig. 13).The mechanical properties of 3DPv are clearly different from those of the 3D printed non-porous particles 3DPf, even though they have the same particle shape characteristics based on the results of image analysis (Fig. 7).Therefore, the experimental results indicate that the 3D printing technique has the potential to control the inner structure in addition to the particle shape.However, few characteristics, such as the OL (Fig. 2b), pore size distribution, anisotropic behavior, and time-dependent behavior, have not been investigated in detail in this study.
In addition to these issues, a notable discrepancy of triaxial compression behavior between 3DPv and VS was observed and needs further discussion.The OLs and pore sizes were randomly allocated because of the input parameters to generate the particles (Table 1 and Fig. 2).These parameters may reflect natural porous particles such as VS, which has an inhomogeneous porous inner structure.However, it appears to lead to undesired variation in the testing results (Fig. 9), even though the variation is considerably smaller than that in VS.According to Ge et al. [10], the modeled aggregate with a consistent BL (Fig. 2a) can result in higher repeatability especially until the first breakage point.Moreover, the widely distributed OL and pore sizes are considered to cause anisotropic behavior of the particle strength under different boundary conditions.Furthermore, it is reported that the anisotropic behavior can be induced from the printing direction to the loading direction [3,10].Therefore, in the study, the experimental results (Figs. 9,10,11,12,13) must be understood as a combined effect of these factors: inhomogeneity of bond strength for smaller unit particles and anisotropic behavior due to the printing method.
In addition, there is a difference in the time-dependent behavior between the 3D printed particles and natural soil.The 3D printed particles exhibit more ductile behavior compared to VS during isotropic consolidation, and the tendency is more significant for porous 3DPv than for 3DPs and 3DPf (Fig. 12).Such behavior of the 3D printed particles has already been pointed out by few previous studies [3].Moreover, the behavior may depend on the 3D printing material [47].More experimental data are required to analyze this issue.
Aside from the issue of time-dependent behavior, a notable discrepancy of volumetric strain behavior in the triaxial test was observed.The amount of volumetric strain under isotropic consolidation for 3DPv is significantly higher than that for VS (Fig. 12), whereas those under shear loading have the opposite tendency; conversely, that for VS is relatively higher than that for 3DPv (Fig. 13).
The results indicate that the 3D printed spherical porous particles (i.e., 3DPv) is more deformable under isotropic consolidation.However, natural porous volcanic soils tend to be crushed more under shear loading rather than under isotropic consolidation [18].The discrepancy may be attributed to the difference of particle shape, fabric anisotropy induced by irregular particle shape and sedimentation process under gravity [32,34], and their relationships with the boundary effects such as loading direction and stress path [46].Further experiment or simulation for porous particles with irregular particle shape will give us a deeper insight into this issue.
However, the use of 3D printed porous particles is expected to enable a better understanding of the mechanical behavior of natural crushable soils such as VS.In this study, the 3D printed porous particles were modeled by an aggregate bonding of the smaller unit spheres (Fig. 1).The porosity results from the volume fraction and diameter of the unit spheres (Table 1).Therefore, these parameters may help in controlling the porosity and maintaining the same particle shape.This may help us reveal the individual effects of porosity on the mechanical properties in a more systematic manner.Another prospect of applying the Fig. 15 Comparison of strength parameter / 0 for various volcanic soils (our results and results of Kawamura [24] for Ta-d (original) in addition to those of Miura et al. [33] for coarse-grained volcanic soils sampled in Mori, Tomikawa, Kashiwabara, Bibi, and Nakashibetsu) proposed method is to model irregularly shaped particles with inner void structures, as is performed in the DEM simulations.Once an irregularly shaped particle is transformed to the STL model, for example, a desired inner structure may be modeled by the aggregate.Moreover, as the 3D printing technique has advanced, the available 3D printing materials have become more diversified, including polymer, gypsum, and sand (with binders).Xu et al. [47] showed differences in the mechanical behavior using gravelly particles printed by polymer and gypsum.Thus, the 3D printing technique enables the creation of controllable experimental parameters related to not only particle shape but also particle strength.

Conclusions
In this study, porous and non-porous particles with identical particle shape were replicated using the 3D printing technique, and their geotechnical properties, such as particle shape characteristics, single particle crushing strength, and triaxial compression behavior, were investigated.The particles with or without inner void structures, namely, 3DPv and 3DPf, were modeled by the DEM simulation program PFC 5.0, which were 3D printed using a resinbased material.For comparison, the 3D printed sphere (3DPs) and natural volcanic soil composed of porous particles (VS) were also tested.The following conclusions were drawn on the basis of the results.
The circumcircle, sphericity, roughness, and roundness were evaluated using image analysis for each material type.The results indicate that 3DPf and 3DPv have the same particle shape characteristics, which are different from those of VS.Because the inner structures of 3DPf and 3DPv are undoubtedly different based on the measured mass of each material, porous and non-porous particles with identical particle shape characteristics can be successfully generated by the proposed modeling method and 3D printing technique.The results of the single particle crushing test proved that the particle strength of 3DPv was significantly lower than those of 3DPs and 3DPf.The difference in inner structures resulted in the difference in the particle strength.Moreover, compared to VS, the variation of the results was significantly small for all types of 3D printed particles.The triaxial compression test results of 3DPv and VS exhibited highly contractive volume changes and stress-hardening behaviors, whereas those of 3DPs and 3DPf exhibited dilative and stresssoftening behaviors.Based on the experimental results, it can be concluded that 3DPv has the typical mechanical properties of porous particles.

Fig. 2
Bond models for aggregates in a Ge et al.[10] and b this study and the c frequency distribution of OL

Figure 4
displays a photograph and X-ray computed tomography (CT) and scanning electron microscopy (SEM) images of the VS used in this study.Figures4a and 4cshow the photograph of the original Ta-d and VS to highlight the physical appearances.The coarsegrained VS with several voids inside the particles (Figs.4b

Fig. 3
Fig. 3 Manufactured 3D printed particle models and particles: STL models of a 3DPs, b 3DPf, and c 3DPv and printed d 3DPs, e 3DPf, and f 3DPv

Fig. 4
Fig. 4 VS particles tested in this study: a original Ta-d; b X-ray CT image of the sampled Ta-d (Diameter of thin-walled sampler = 75 mm).The black parts inside particles are inner voids.c tested material (after oven-dried and sieved to 2-4.75 mm), and d SEM image for a particle of VS

Fig. 6
Fig. 6 Schematic of the single particle crushing test apparatus Figure 9a and b indicate that 3DPs and 3DPf typically follow three distinct axial force vs. axial displacement responses with different yield strengths and corresponding

Fig. 9
Fig. 9 Axial force and axial displacement responses of a single particle crushing test of a 3DPs, b 3DPf, and c 3DPv

Fig. 11
Fig. 11 Failure mode on single particle compression of a 3DPs (Left: before compression, right: 0.41 mm in axial displacement), b 3DPf (Left: before compression, right: 0.75 mm in axial displacement) c 3DPv (Left: before loading, right: 0.53 mm in axial displacement), d VS that exhibited compression failure (Left: before compression, center: 0.24 mm in axial compression, right: after compression), and e VS that exhibited tensile failure (Left: before compression, center: 0.33 mm in axial compression, right: after compression)

Table 1
List of parameters for particle generation

Table 2
List of parameters for the 3D printed particle models

Table 3
Properties of STL

Table 5
Physical and mechanical properties of the original Ta-d sampled in Hokkaido Physical properties Values Natural water content, % 144.6-190.9b Tested based on JIS A1202, Test method for density of soil particles c Based on the model size (Also see Tables 1 or 2) d

Table 7
Testing conditions of triaxial specimens