Experimental and numerical determination of representative elementary volume for granular plant materials

A series of physical and numerical tests were conducted to determine representative elementary volume of granular plant material. The load response of pea grain assembly poured into a cuboid test chamber and subjected to uniaxial confined compression was studied. The apparatus was equipped with adjustable side walls that allowed measurement of boundary stresses in samples of varying thickness. It was found that load distribution varied considerably in samples of thickness smaller than three times the size of the particle. Less pressure variation was observed in grain assemblies of thickness equaled to three, five and seven times the particle size. Comparison between experimental data and numerical DEM results have shown qualitative agreement. It was found that the specimen of dimension not smaller than five times the particle size can be used as a representative elementary volume in confined uniaxial compression test of granular plant materials.


Introduction
The large-scale application and processing of granular plant materials in agriculture and many branches of industry (e.g. cosmetic, pharmaceutical or food industry) requires improved insight into the complex nature of grain and seed assemblies. Experimental and analytical research methods have been used for many years to investigate and predict phenomena occurring in grain assemblies. The difficulties in the interpretation of many processes observed in granular solid resulted in development of computational techniques which are useful for studying the behaviour of the grains. Both physical and numerical tests have their disadvantages. Full-scale physical tests are laborious and expensive while the numerical ones require time-consuming calculations and high computational power, especially when computational techniques based on micromechanical approach are used. One of the most popular methods allowing the discrete nature of granular material to be taken into account is the discrete element method (DEM) [6].
One of the important questions regarding examinations of granular materials concerns the representative elementary volume (REV) which is the smallest specimen over which a measurement can be made that will yield a value representative of a large volume (bedding). The choice of representative elementary volume is of high importance in the study of materials. Even in the case of homogenous solids the answer is not obvious. Gitman et al. [11] proposed that for elastic materials the size of the REV can be determined, however for materials showing localized deformation in the regime of softening, the REV may not exist. The size of representative volume of granular matter can be associated to the size of grain and it may be determined analytically, numerically or based on experiments. The size of specimen has been shown to have a significant influence on the mechanical response of granular solid. Based on extensive past experience, the authors of Eurocode 1 [9] recommend the diameter for the direct shear cell to be 20 times of the largest particle size. Size of sample affects load transmission and fabric in grain assembly [27,31]. Masson and Martinez [18] showed that sample size equal to 7 to 8 times the size of the largest particle was sufficient to obtain a macroscopically representative value when porosity and coordination number were analyzed whilst a size 12 times the largest particle size was required for proper stress tensor evaluation. These results were in agreement with experimental data of Lanier et al. [16] and Biarez and Hicher [3].
A comparison of the load-displacement responses obtained for rapeseed samples of various thickness under uniaxial compression showed some irregularities in the load-displacement relationships for specimens of thickness lower than 11 times the largest particle size [27]. These authors reported that an increase in the specimen size resulted in an increase in the observed lateral-to-vertical pressure ratio.
The difficulty with the choice of REV in numerical simulations can be reduced through application of periodic boundaries which allows for modeling sample free from boundary effects [17,22,26,30]. When the center of particle protrudes outside the periodic boundary, a mirror particle is generated on the opposite side of sample to maintain continuity in that dimension. Martin et al. [17] simulated uniaxial compaction of bronze, copper and aluminum powders composed of spherical particles. Simulations were run for periodic cells composed of 400 or 4,000 uniform particles. The difference between the two samples came out only in the smoothness of the resulting curves.
The brief review of previous studies has shown that there is little consensus on how sample size affects the observed or predicted properties of a granular solid. In this paper, the influence of sample size on the mechanical response of uniaxially compressed granular assembly was studied to establish the representative elementary volume in this loading situation.

Uniaxial confined compression test
Uniaxial confined compression test is a common method to determine mechanical properties of granular materials which are of interest to technological process designers. Both physical testing [27,29] and numerical [1,30] modeling of compression of granular solid can provide valuable knowledge necessary for efficient design and useful scientific insight.
Uniaxial compression test allows one to analyze the stressstrain characteristics of granular solid and to determine material parameters such as lateral-to-vertical pressure ratio k, modulus of elasticity E or Poisson's ratio ν. Eurocode 1 [9] recommends the measurement of the lateral pressure ratio k, defined as the ratio of horizontal stress to the vertical stress, in specimen of diameter D larger than 5 times the largest particle size and height in the range from 0.3D to 0.4D. The recommendation of Eurocode 1 enhances the repeatability of test. The Eurocode, as well as the majority of published analyses considered cylindrical sample under axi-symmetrical stress with the assumption of material isotropy. Numerous examinations have pointed to considerable anisotropy of granular materials dependent on the shapes of particles and method of generation of the bedding [20,24]. Therefore for this study the uniaxial apparatus of rectangular cross section was constructed to measure loads in the two perpendicular horizontal directions x and y in the sample undergoing deformation in the vertical direction z. The construction of apparatus was far from the one recommended by Eurocode 1, however it should not have a significant impact on the uniaxial compression test. The granular material poured into the test chamber was compressed through the top platen that moved vertically with a constant velocity.
The curve illustrating typical load-displacement relationship includes two phases: loading (A-B) and unloading (B-D) (see Fig. 5). The nonlinear loading path includes both elastic and plastic deformation of specimen. The initial phase of unloading (B-C) is nearly linear and a predominantly elastic response, followed by nonlinear path where elastic and plastic strains take place in the sample simultaneously [29].
It is generally accepted that the pressure ratio is influenced by the microstructural properties, the packing density and resulting fabric of the granular solid. Horabik and Rusinek [13] showed the influence of moisture content and grain shape on the stress distribution in a sample subjected to uniaxial compression. They reported the pressure ratio decreased with increasing sphericity of grains. The k value was measured to be 0.37 for nearly spherical soya beans, while k of 0.74 was found for disk shaped lentil grains. An increase in moisture content of grains was also found to result in a decrease in pressure ratio.
A study of the heterogeneity of load transmission in rapeseed samples compressed in a cuboid chamber, conducted by Rusinek et al. [27], revealed strong relation irregularities in the load-displacement response as the thickness of sample was varied. Comparison between DEM predictions and physical load-displacement responses of glass beads and corn grains subjected to compression test by Chung et al. [4] showed that DEM results were in good agreement with experimental data only for non-spherical grains. DEM appeared to overpredict the forces on the bottom platen for the spherical glass beads, probably as an effect of representation of real glass beads using perfect spheres in numerical simulations.
The laboratory tests were conducted in this study using a uniaxial compression apparatus of rectangular cross-section The confined uniaxial compression tests were conducted for chicken pea (Piast variety) of geometric mean diameter of 7.4 mm [19] at moisture content of 11 %. Three axial dimensions of 50 seeds were measured with major a, intermediate b and minor c dimensions having mean and standard deviation as 8.00 ± 0.30, 7.30 ± 0.26 and 6.90 ± 0.38 mm. Pea was taken as experimental material because of its nearly spherical shape and dimensions allowed for reasonable design of experiment and DEM simulation. The seeds were poured into chamber with extremely eccentric stream against a side wall. The surface of specimen was then leveled before starting the compression test.

Discrete element method
The discrete element method (DEM) is a common numerical technique for detailed investigation of granular systems based on a microstructural approach [6]. The method, originally proposed to model interactions between rigid rock blocks and to study ground mechanics is nowadays a promising tool for modeling 2D and 3D processes in granular solids [5,10,12].
Detection of contacts between particles is followed by calculating contact forces at each incremental time step. The time step is set small enough to allow one to assume a constancy of translational and rotational accelerations. Integration of equations of motion for each particle provides velocities and positions of particles. The rigid particles are allowed to overlap locally at contact points using a soft contact approach.
In this study a simplified non-linear Hertz-Mindlin contact model [14] with elastic spring and viscous damper in the normal direction (see Fig. 2a) and spring, damper and frictional  slider in the tangential direction (see Fig. 2b) was applied. Spring models accumulation of elastic energy in the system, whilst damper and slider model the energy dissipation. The tangential contact force is limited by the Coulomb friction law which assumes that particles slide over each other when the tangential force is at limiting friction. Discrete element method allows for detailed study of, inter alia, compression [1,4], shearing [12,28] and milling [23] of granular materials. In this study, three-dimensional DEM simulations were conducted using the EDEM software. Uniform 7.3 mm diameter spheres (d) with mechanical parameters of pea grains were poured into a steel chamber of rectangular cross-section and various thickness (T ) (see Fig. 3). Experimentally established input parameters for the grains [31] and the steel [32] are listed in Table 1.
The modelled assemblies were generated by filling the chamber with 260, 770, 1,200 and 1,630 monosized particles. Particles randomly generated in the whole volume of box settled down onto the bottom of box under gravitational forces. Initial configuration of 21 mm thick numerical specimen is presented in Fig. 3. The top lid of box was then moved down with constant velocity of 0.3 × 10 4 mm min −1 until the normal lid pressure reached 100 kPa. The unloading was modeled by moving the lid up until there was no contact between specimen and platen. The factor multiplier of 0.1 (10 %) was applied to the computed critical time step for DEM simulations, ensuring numerical stability without increasing computational cost [25]. Three replications were performed for each thickness of the chamber.

Numerical experiment
The effective elastic modulus E was used to characterize the bedding elasticity [8] and was calculated for each 1 μs time step as follows: where H is the height of the assembly, Δσ z is the change in vertical pressure during each time step, ν is change in vertical displacement and k L is a ratio of the change in lateral pressure to the change in vertical pressure during each time step.
The porosity for 1 diameter thick sample was found to be 0.53 which is distinctly higher than those for thicker samples all at approximately 0.48. The modulus E was found to increase with σ z increasing from 0 to 100 kPa ( Fig. 4) with some fluctuations that were larger for thinner samples. No significant difference in the overall E(σ z ) characteristics for various thicknesses of the sample was found. However an increase in thickness T resulted in a decrease in the standard deviation due to the lower scatter resulting from averaging over a higher number of particles.
The comparison between the lateral-to-vertical pressure ratio in the x and y directions (see Fig. 1) for specimens of various thicknesses showed that there was no significant difference between the pressure ratios calculated for pressures directed along the x axis (k x ) (see Fig. 6a). The change in thickness of the sample in the y direction from 1 to 7 seed diameters did not influence significantly the pressure acting in the x direction. Only the curve for one seed diameter thick sample was located slightly above the three remaining probably due to the weakest averaging of contact forces acting in a small number of contact points. In the case of 1 diameter thick sample much fewer contact points existed at the 0.12 m height of the sample causing very uneven and unrepeatable distribution of contact points between seeds and with the walls. Contrary to the pressures in the x direction, the pressures calculated in the y direction varied strongly with sample thickness and so influenced k y accordingly (see Fig. 6b). Each step of increase in sample thickness from 1 to 5 seed diameters resulted in an increase in the number of contact points. As a result approximately sixfold increase in k y (from 0.074 to 0.473 at 50 kPa of vertical pressure) was observed when the sample thickness increased from 1 to 5 seed diameters. Further increase of thickness from 5 to 7 seed diameters did not bring noticeable increase in k y . Following this result numerical sample width equal to five particle diameters may be recommended as minimum sample width to obtain repeatable values of pressure ratio. The comparison between the lateral-to-vertical pressure ratio in the x and y directions showed that k y values was comparable to k x values above a thickness larger than five seed diameters.

Comparison between numerical and experimental results
The DEM predictions are compared with experimental results in Figs. 4, 5 and 6. The results are plotted using the mean values with the error bars indicating ± one standard deviation. The initial porosities (Φ), effective elastic moduli (E), the load-displacement responses and lateral-to-vertical pressure ratios of numerical and physical samples under compression load were compared. The porosities for numerical samples thicker than one particle diameter were found to be all approximately equal to 0.48 while those of the larger physical specimens composed of non perfect spheres (pea grains) were all approximately 0.53 or 10 % higher. The uniform spheres used in the DEM simulations gave a closer packing and a lower porosity [12].
As shown in Fig. 4, the DEM simulations predicted a lower effective elastic modulus than the experiments for the thinnest sample. Two factors might have major influence on the observed difference in material behavior in the thinnest layer. One is the not perfectly spherical shape of pea seeds and the other is the crossover from quasi-2D behavior to 3D behavior. Both of these factors result in increase in the strength of the bedding. Results of Azema et al. [2] and of Donev et al. [7] gave some explanation of the role of particle shape. Azema et al. [2] observed that in granular assemblies composed of nearly circular particles, slight deviations from disc shape have strong space filling effect. The authors found that in the 2D assemblies the shear strength was a linear function of particle elongation. Increasing mobilization of friction and associated anisotropy were proposed to be key effects arising from particle elongation. From 3D simulations with spheroids, Donev et al. [7] reported that the number of contacts was the lowest for perfect sphere and rose sharply for small deviations. They found that the maximum coordination number was reached at particle aspect ratio λ of approximately 1.3. Landry et al. [15] conducted large scale DEM simulations of granular packings in 2D, quasi-2D and 3D to study the crossover from two to three dimensions. In simulations of packings of spherical particles in rectangular containers a width of two sphere diameters was found to represent a transition to more 3D-like behavior. The authors focused on the vertical stress in simulated packings and observed that 2D packings supported much less vertical stress than 3D packings. It was also shown that polydispersity and different force models did not have a large effect on the vertical stress in 3D packings.
The difference between stiffness of numerical and physical samples decreased significantly for specimens of thickness from 3 to 7 times particle size. For these samples, the difference between numerical and experimental E values lie within the range of scatter which is probably a result of combination of the imperfections of real granular material and the imperfections of experimental apparatus as compared to configuration of numerical experiment. Imperfections of real apparatus include among others: limited accuracy of its dimensions, non perfectly rigid walls and their supports, imperfect plane surfaces or uneven distribution of frictional properties. These effects had less influence in larger samples. The behavior of larger samples of real material thus became closer to that of ideal material used in the numerical experiments.
The study has shown that the mechanical response of granular assembly subjected to uniaxial compression is significantly affected by dimensions of the system. Both numerical and physical tests revealed difference in stress transmission in spherical or nearly spherical assemblies of thickness lower than five times the size of particle. The findings suggest that a specimen of dimension not lower than five times the size of the particle should be used as a representative elementary volume (REV) for uniaxial compression testing. Figure 5 shows the relationship of vertical pressure versus volumetric strain for the numerical and experimental samples of 33 mm thick (five times particle size). The both, numerical and experimental curves illustrate that the average trends of three replications are qualitatively similar. Initial slopes in the plots are approximately equal, but beyond ε of approximately 0.008, inclination of the DEM curve became steeper.
Numerical sample was less compressible. It reached 100 kPa of vertical stress at volumetric strain of 0.019 while experimental sample reached the same value of σ z at ε of 0.021 or approximately 1.1 times higher. The first phase of unloading was found to be approximately linear and of nearly the same inclination in the two relationships. Such a linear elastic relationship during unloading of sample of grain was reported by Sawicki [29] and by Rusinek et al. [27]. The last part of unloading curve was steeper in the case of DEM sample. Values of residual volumetric strain after unloading was 0.0093 in numerical and 0.012 in experimental samples respectively. Eurocode 1 [9] recommends using the linear part of the pressure-displacement response during unloading to calculate the effective elastic modulus. The effective elastic moduli during unloading were calculated as 4.27 and 4.54 MPa for the experimental and numerical samples respectively. DEM overpredicted slightly the stiffness of sample as an effect of representation of real grains using perfect spheres in numerical simulations and application of a contact model describing real contact only approximately. It is commonly known that the mechanical parameters of plant materials measured experimentally can vary widely due to the lack of standardization in determining their values and the variation of properties of materials caused by conditions of cultivation and plant variety. One such examination showing wide variability of mechanical parameters of plant materials measured in two distant laboratories in Spain and in Poland was that of Molenda et al. [21] where effective elastic modulus of uniaxially compressed wheat varied from 5.79 to 22.4 MPa.
The comparison between numerical (Fig. 6b) and experimental (Fig. 6c) lateral-to-vertical pressure ratios in the y direction show that DEM predicted two times higher pressure ratios for the thinnest samples. The pressure ratios of 0.07 and 0.14 were calculated from physical and numerical data, respectively, at vertical pressure of 50 kPa. This difference is due to the two factors: not perfectly spherical shape of pea seeds and the crossover from quasi-2D behavior to 3D behavior, which was described at Fig. 4a. The pressure ratios of 0.32 were calculated for samples of 21 mm thick. For the samples of 33 and 44 mm thick the differences between numerical and experimental k y value lie within the range of scatter.

Conclusions
The influence of specimen size on the mechanical response of granular assembly subjected to uniaxial compression has been investigated using DEM and physical experiment. The study was conducted using pea grains placed in a chamber of rectangular cross-section.
The chief conclusions are as follows: • The effective elastic moduli of both numerical and experimental assemblies were not affected by sample thickness of 3, 5 and 7 diameters of particle. • The lateral force transmission in confined compression in the form of ratio of horizontal to vertical pressure was found to be sensitive to sample size when the sample thickness is less than five times the particle size. • A good qualitative agreement between the numerical and experimental vertical pressure-volumetric strain relationships and effective elastic moduli was obtained for granular specimens of various sizes. • DEM predicted a slightly stiffer response for the spherical assembly as compared to the experiments on pea grains. The effective elastic modulus was found to be 4.27 and 4.54 MPa for the experimental and numerical samples respectively, for sample thickness of five times particle size. • A specimen of dimension greater than five times particle size can be used as a representative elementary volume (REV) in uniaxial compression confined test of granular plant material.
The above conclusions were drawn based on representing the pea grains using spherical particles. Determining the REV for non-spherical particles requires further investigation.