Packing, alignment and flow of shape-anisotropic grains in a 3D silo experiment

Granular material flowing through bottlenecks like the openings of silos tend to clog and to inhibit further flow. We study this phenomenon in a three-dimensional hopper for spherical and shape-anisotropic particles by means of X-ray tomography. The X-ray tomograms provide information on the bulk of the granular filling, and allows to determine the particle positions and orientations inside the silo. In addition, it allows to calculate local packing densities in different parts of the container. We find that in the flowing zone of the silo particles show a preferred orientation and thereby a higher order. Similarly to simple shear flows, the average orientation of the particles is not parallel to the streamlines but encloses a certain angle with it. In most parts of the hopper, the angular distribution of the particles did not reach the one corresponding to stationary shear flow, thus the average orientation angle in the hopper deviates more from the streamlines than in stationary shear flows. In the flowing parts of the silo shear induced dilation is observed, which is more pronounced for elongated grains than for nearly spherical particles. The clogged state is characterized by a dome, i. e. the geometry of the layer of grains blocking the outflow. The shape of the dome depends on the particle shape.


Introduction
Silo discharge and the flow of granulates through narrow orifices are important ubiquitously in everyday life and in technological processes. A huge number of experimental, as well as numerical studies have been performed over many years to understand the underlying processes (e. g. ).
The spontaneous formation of arches [19,20], the kinetics preceding arch formation [22], forces and force chains [23,24], and arch breaking [11][12][13] are among the important features when one is interested in the prevention of silo clogging. The formation of arches, i. e. arrangements of grains at the outlet of a two-dimensional hopper that block outflow, was analyzed in a pioneering study by To et al. [14]. For disks in a two dimensional silo, they found that the arch span is always larger than the opening of the silo and that for frictionless particles the arch is convex everywhere. Introducing friction they found clogging events where the arch at the opening may not be convex locally. In addition, they determined the jamming probability as a function of the size of the hopper opening. Clogging probabilities were also measured not only for spheres but also for one type of elongated grains by Zuriguel et al. [25] in a three dimensional silo. In their model, the mean avalanche size between two clogging events could be fitted with a power law. The model assumes a critical radius R c of the orifice, above which the hopper outflow is continuous, without clogging. However, the divergence of the mean avalanche size with increasing orifice radius R according to the model in Ref. [25] is proportional to (R c − R) −7 , so it is very difficult to confirm the existence of such a critical radius experimentally with reasonable amounts of grains available.
In the present experimental study, we analyze the effect of shape anisometry of granular particles on hopper flow and clogging by means X-ray tomography. For this purpose we use elongated pegs, flat lentils and nearly spherical peas. These particles have sizes of several millimeters, and can be detected and individually distinguished in the recorded X-ray tomograms. The silo is a cylindrical bucket with changeable orifice size. The size is chosen such that clogging occurs after finite avalanche sizes. The tomograms are recorded in the clogged state with the grains at rest. From comparison of subsequent clogs, we can characterize the flow field (streamlines) of the granulate in the silo. We find all positions and orientations of individual particles in the clogged state, and use the tomographic data to calculate the packing fractions and particle orientations. The tomograms also give information about the shape of the domes, as well as on the orientation of particles forming the dome.

Experiment
The experiments were performed with a small plastic silo consisting of a bucket of a diameter of 19 cm and a height of 21.4 cm. The aperture of the outlet was adjustable by exchanging insets in the center of the bottom plane. The bottom of the bucket has a quadratic opening of about 5 × 5 cm 2 , where different plates with circular holes were inserted, hole sizes could be chosen in 0.5 mm radius steps. A schematic drawing of the silo is shown in Fig. 1. The different granulates studied in this setup are also shown in Fig. 1. To characterize the deviation from the spherical shape we use two different quantities: (1) the aspect ratio Q, which is the ratio of the size of the particle along its rotational axis and perpendicular to it and (2) the equivalent radius r eq which is the radius of a sphere with the same volume as the particle. The following materials were used: • wooden pegs with cylinder shape and tapered ends, diameter d = 5 mm, length = 40 mm, aspect ratio Q = /d = 8, equivalent radius r eq = 5.7 mm • wooden pegs with cylinder shape and tapered ends, d = 5 mm, = 25 mm, Q = 5, r eq = 4.9 mm • wooden pegs with cylinder shape and tapered ends, d = 6 mm, = 20 mm, Q ≈ 3.3, r eq = 5.1 mm • wooden pegs with cylinder shape, d = 5 mm, = 10 mm, Q = 2, r eq = 3.6 mm For the preparation of the experiment, the silo outlet is closed first. Then the container is filled with the grains. We have selected the outlet sizes for the individual grain types such that after the bottom hole is opened, the granulate flows out but clogs after an appropriate time. The conditions for a well prepared state are that the avalanche is sufficiently large to create a representative flow and alignment pattern in the container, and that the avalanche stops by clogging when there is still enough material in the container cover more than the region covered in the subsequent tomogram. In this clogged state, an X-ray tomogram of the lower half of the silo is recorded. We use the robot-based flat panel X-ray C-arm system Siemens Artis zeego of the STIMULATE-lab, Otto von Guericke University, Magdeburg. The chosen spatial resolution was 2.03 pixel/mm, with recorded volumes of 25.2 cm × 25.2 cm × 19 cm. The 3D arrangement of the particles is then determined from the X-ray tomogram.
To calculate the packing density in different parts of the silo, the tomograms are binarized. We define a threshold for all voxels, that determines whether the voxel belongs to a grain or not. For that we have to choose carefully the threshold. Our algorithm is based on Otsu's method [44] and to find a good threshold, we only take the part of the tomogram with particles in it. Borders and empty regions are avoided. Afterwards the processed tomograms are averaged for each material. After this the 3D image is projected to a 2D representation where we average over all equivalent voxels, i. e. voxels with the same distance to the center of rotation and the same height. This method gives a measure of the relative packing densities in different parts of the silo.
In order to calculate the streamlines, the nematic order parameter and the average orientation angle, we determined the position and the orientation of the individual particles; for details see [45]. The shearinduced orientational order is monitored by diagonalizing the symmetric traceless nematic order tensor T : where (n) is the unit vector along the long axis of particle n, and the sum is over all N detected particles in equivalent regions of the silo. For that purpose, summation over a large number of particles is needed. As a compromise, we make the assumption that the flow in the axially symmetric geometry of our silo yields an axially symmetric mean flow field and consequently an axially symmetric mean arrangement of the particles, therefore, we average over ring-shaped zones of the silo. Conventionally, the nematic order parameter S is the eigenvalue of the order tensor with the largest absolute value. In the case S > 0 (preferential alignment along a certain axis), we will refer to S in the definition of flow alignment and define the average flow orientation angle θ av as the angle betweenê S (the eigenvector of T corresponding to S) and the local flow direction, see Fig. 2

(a).
In regions with strong planar influences, i.e. near walls and close to the free surface, the orientational ordering of the particles is often fan-like: the axes of the particles are preferentially aligned perpendicular to a certain axis, but azimuthally more or less distributed at random. In that case, T has two small positive eigenvalues, and the third one (S) is negative and has the largest absolute value (recall that T is traceless). In that case, we use a different, more descriptive concept to characterize the local alignment: The eigenvectorê S corresponding to S marks the direction in which the particles are preferentially not oriented. Therefore, θ av is defined differently, as the angle between the streamlines and the closest vector in the orientational plane (plane perpendicular toê S ), see Fig. 2(b). In this case, we characterize the flowinduced order by the largest positive eigenvalue S . This helps us to visualize shear induced orientational changes in the entire hopper, since S is less susceptible to boundary effects.

Order, Orientation and Packing Fraction for Pegs
Central vertical cuts of two subsequent tomograms are presented for pegs with aspect ratio Q = 5 in Fig. 3. These tomograms were taken after several avalanches, thus the orientation and packing of the  grains reflects the properties of the developed flow (initial conditions erased in the flowing regions). As explained above, these tomograms were taken in the clogged state. During the avalanche between these two tomograms, approximately 60 particles left the hopper. The displacement of the material due to this avalanche is visualized by a few selected particles which are marked with colors in Fig. 3. It is clearly visible that the central part moves faster, especially approaching the outlet. Stagnant zones can be identified near the walls especially at the bottom of the hopper. In the flowing region particles are expected to be oriented nearly (but not exactly) parallel to the streamlines [37,46,47]. In order to quantify the ordering and alignment of the particles and determine the flow field and the density distribution, tomograms from 26 subsequent clogged states have been analyzed. In Figures  4a-c, the colormaps indicate the density distribution, the effective order parameter (S ) and the average alignment angle (θ av ) of the particles in the r − z plane. The alignment angle is defined in the interval −90 • < θ < 90 • . The alignment angle is measured with respect to the streamlines, which are indicated by white lines. In order to better visualize the correlation between these quantities, on Fig. 4a we represented the two other parameters (the average orientation of the particles and the nematic order parameter) by ellipses. The size and orientation of the ellipse represent S and θ av in the following way. The area of the ellipse is proportional to the projection of the average orientation of the particles to the r − z plane. In the flowing region the rods are aligned due to shear in this plane, so we find ellipses with large area. Near the vertical walls, however, the particles are mostly oriented tangentially, thus the area of the ellipse is small. The direction of alignment is reflected in the orientation of the ellipses. The order parameter is represented by the flattening of the ellipses: S = flattening = 1 − b/a, where "a" and "b" are the lengths of the long and short semi-axes, respectively. Thus S = 1 (perfect order) would be represented by a thin line, and S = 0 (isotropic state) by a circle. All the above applies for usual (calamitic) nematic order, S > 0, where the largest absolute-value eigenvalue of T coincides with the largest eigenvalue (full ellipses). For the fan-shape oriented case, where the eigenvalue with largest absolute-value is negative (S < 0), we draw an empty ellipse with flattening |S| = 1 − b/a. In those cases the short (of length "b") semi-axis is oriented to the axis of the corresponding eigenvector of the order tensor. At these locations, the history of the filling procedure or boundary effects produced an initial 2D order, where the orientation of the particles are nearly within a plane, resulting in a negative order parameter. To quantify the flow induced ordering it is more suitable to take the largest positive eigenvalue as an effective order parameter S .
Altogether, Figure 4 shows that in the flowing (sheared) regions we find well oriented rods due to the shear flow, while near the walls tangential alignment is preferred due to the boundary conditions. We see, that in the shear flow the average orientation is not parallel to the streamlines, but there is a flow alignment angle as expected [37,46,47]. The sheared regions are characterized with lower packing fraction, consistent with our earlier observation in shear flow, where the shear induced dilation was found to be only partially compensated with a slight density increase due to the ordering of the grains [46].
In order to accurately characterize the alignment of the particles in this complex flow field, let us recall earlier data on the behaviour of these particles in a split bottom Couette flow [48]. In that experiment, we sheared the same material in a cylindrical geometry and recorded how the orientational order develops from an initially random sample. This dataset is presented in Fig. 5a by a black line. As the local shear deformation γ is increasing, the order parameter rapidly increases and by γ = 1, it already reaches about S = 0.6. During this time, the average alignment angle strongly decreases from about θ av = 80 • to about θ av = 35 • as a consequence of the continuous rotation of the rods due to the shear flow. When the sample is sheared further (above γ = 6), a stationary state is reached which is characterized by a shear alignment angle of about θ av = 8 • . In this state the particles still rotate, but their rotational velocity strongly depends on their orientation, with slow rotation around θ = 8 • , resulting in an ensemble averaged alignment of the sample. This dataset was obtained from 90 X-ray tomograms in the cylindrical split bottom geometry.
Focusing on the orientation of the particles in the hopper (see colored data points in Fig. 5a), we find that a good part of the sample is well ordered with similar values of S and θ av (see blue data points) as in stationary shear flow. This is the region right above the outlet (indicated with blue) in Fig. 5b, where the sample was exposed to the largest shear deformation. Further upstream where the material was sheared less, we find a smaller order parameter and a larger alignment angle, similarly to the initial stages of the case of simple shear (see regions with light blue, green, yellow and red colors in Fig. 5a-b). In order to better track the evolution of the particle orientation in the hopper flow, we present the histogram of the particle orientations for 3 regions in our hopper in Fig. 5c. The three regions are marked with different colors (see the inset of Fig. 5c). The corresponding 3 curves show gradual narrowing of the angular distribution and decreasing of the average angle as the material is exposed to larger and larger deformation as it is displaced downwards in the hopper. In Fig. 5c we also show the histogram of angles of stationary split bottom Couette flow (black line). We see that near the outlet (blue curve) the angular distribution almost reached the case of stationary shear. Coming back to Fig. 5a, it is remarkable, that the data points taken in the hopper are almost all below the black line (corresponding to the case of simple shear). Physically, the difference between the two cases is, that in the hopper the streamlines are strongly converging while in the cylindrical split bottom geometry they run parallel with each other. Thus in the hopper the converging flow leads to particle orientations which are closer to the streamlines than in a simple shear flow.

Packing Fraction for Various Materials
So far we have presented the flow and orientation fields and the packing density of the particles for one type of pegs (Q = 5) in a hopper. For this we used 26 tomograms, with small particle displacement (small avalanches) between the tomograms, so that the flow lines could be identified. In the following, we present the packing density distribution in a hopper for six materials. For these plots, a similar number of tomograms was taken for each material, but here we allowed more grains to flow out of the hopper between subsequent measurements (compare outlet diameters) to reduce correlation in the arrangement of the grains. The resulting packing density plots are shown in Fig. 6 for pegs with aspect ratio Q = 2, Q = 3.3, Q = 5 and Q = 8, for lentils and peas.
As it is seen in Fig. 6, the density is reduced in the sheared region (compared to the density in the stagnant zone near the walls) for all materials due to the shear induced dilation. As expected, the density in the stagnant zone φ decreases with increasing aspect ratio for all rodlike grains. The region where φ was determined is marked with a black rectangle on the images. The shear induced dilation is stronger for rods, compared to the case of nearly spherical particles (peas). Approaching the outlet from the top along the symmetry axis of the hopper we find a decreasing packing density. Near the outlet there is a practically empty region (marked with blue), the region below the dome. The upper border of this region corresponds to the dome which holds the material above in this clogged configuration. We indicated this with a black line, corresponding to 50% of the packing fraction observed in the stagnant zone. The contour of a hemisphere of equal diameter with the outlet is also marked (as a white line) to visualize deviations from a spherical shape. Analyzing the shape of the dome, we find that its height decreases with particle elongation for the sequence Q = 2 → 3.3 → 5. For lentils and nearly spherical peas the dome is similar to the case of Q = 3.3.
For pegs with Q = 8 we find a higher dome, we note however, that these particles are relatively large compared to the system size, thus boundary effects might already have an influence in this case. We should also note that boundaries also induce layering of the particles especially near the bottom wall, where anisometric particles lay flat. In the stagnant zones near the boundaries the inhomogeneities of the density field are stronger compared to the flowing regions, since there is less averaging over different configurations.

Summary
We presented experimental data about the packing fraction, grain alignment, orientational order parameter, and flow field in a 3D hopper based on X-ray CT measurements. We analyzed subsequent clogged states for 6 materials including elongated particles (pegs), lentils, and nearly spherical grains (peas). We have shown that for elongated particles the grains get ordered in the flowing parts of the silo. Similarly to the case of simple shear flows the average orientation of the rods is not parallel to the streamlines but encloses a small angle with it. The order parameter increases as the grains travel downwards the silo and the local shear deformation grows. In most parts of the hopper the orientational distribution of the grains did not reach the stationary orientational distribution observed for simpler shear flows. Consequently, in these regions the order parameter is smaller and the average alignment encloses a larger angle with the streamlines than in a simple shear flow. The packing density of the material is smaller in the flowing (sheared) regions compared to the stagnant zones. This density decrease is stronger for elongated grains compared to the case of nearly spherical particles (peas). Along the vertical symmetry axis of the hopper the density continuously decreases from top to bottom.
Near the outlet the density distribution averaged over many clogged configurations provides information about the shape of a typical dome, i.e. the last layer of the grains which blocked the flow. A characteristic feature is the ratio of hole diameter and elevation of the dome above the outlet in the center. We find that the height of the dome decreases with increasing aspect ratio for pegs with Q = 2 → Q = 3.3 → Q = 5. This has consequences for the average number of grains that are involved in the formation of the clog, and therefore also for the clogging probability, the avalanche statistics and the critical radius for the divergence of the mean avalanche sizes. The analysis of these details requires not only the information on the dome shape, but also on the characteristic alignment of particles forming the clog. This will be treated in a subsequent investigation.