Modelling the wall friction coefficient for a simple shear granular flow in view of the degradation mechanism

Abstract A steady granular flow experiment was performed in a confined annular shear cell to examine how the wall friction coefficient $\mu _w$ degrades from the intrinsic sliding friction coefficient $f$ between the grains and the container wall. Two existing models are invoked to examine the decay trend of $\mu _w/f$ in view of the ratio of shear velocity to the square root of granular temperature $\chi$ (Artoni & Richard, Phys. Rev. Lett., vol. 115, 2015, 158001) and the ratio of grain angular and slip velocities $\varOmega$ (Yang & Huang, Granul. Matt., vol. 18, issue 4, 2016, p. 77), respectively. As both models correlate $\mu _w/f$ to different flow properties, a hidden relation is speculated between $\chi$ and $\varOmega$, or equivalently, between the granular temperature and the grain rotation speed. We used experiment data to confirm and reveal this hidden relation. From there, a unified $\mu _w/f-\chi$ model is proposed with physical meanings for the model coefficients and to show general agreement with the measured trend. Hence we may conclude that both the fluctuations in grain translations and their mean rotation are the crucial yet equivalent mechanisms to degrade $\mu _w/f$.

While the hydrostatic p has been confirmed in both experiments and discrete element (DE) simulations (Artoni et al. 2018;Lin, Jiang & Yang 2020), DE simulation reveals a non-constant μ w that decays from the sphere-wall sliding friction coefficient f towards the creeping regime in steady, simple shear and surface flows (Richard et al. 2008;Brodu, Richard & Delannay 2013;Artoni & Richard 2015b;Artoni et al. 2018).To explain the weakening, Richard et al. (2008) consider how intermittent and random sphere fluctuations render individual sphere-wall friction events in different directions that cancel each other out to give lower bulk friction.Later, Artoni & Richard (2015b) defined χ = u / T , using the shear velocity u and the streamwise bulk granular temperature T , to indicate that a flow at small χ experiences more incoherent grain dynamics.How μ w degrades from f has been fitted, with two constants (Artoni & Richard 2015b), as (1.1) From a different aspect, Shojaaee et al. (2012) simulated two-dimensional shear flows to study how individual disk rotation may destroy pure sliding motion at the contact points to render μ w lower than f .Such a rotation-induced degradation of μ w is also reported from a simulated three-dimensional granular avalanche by Yang & Huang (2016), together with a degradation model developed in view of the single-grain dynamics as follows.When a grain moves against the lateral wall at translation velocity u pi and angular velocity ω pi , it would develop a total velocity at the contact point as u pi + R × ω pi , with R denoting the radius of gyration from the grain centre.A non-zero ω pi would divert the surface sliding velocity, together with the generated Coulomb friction, away from u pi .Hence the freedom of individual grain rotation provides a different friction-cancelling mechanism at the grain-size level that can accumulate to degrade μ w mechanically according to the degree of rotation-induced diversion.The authors proposed a dimensionless rotation index Ω = |ω| D/2u with mean angular speed |ω| = |ω pi | , and sphere diameter D. They proposed a different model as with its coefficients a = ω ⊥ /|ω| and b = ω /|ω| depending on the relative strength of the mean angular velocity components perpendicular and parallel to u , respectively.In fact, Louge (1994) also studies the sphere-wall friction force for a rapid granular flow in simple shear.After normalizing the tangential stress by the normal stress at the wall, a non-constant μ w /f is found to correlate well with a dimensionless relative velocity at the contact plane, r = (|u | + |ω| D/2)/ 2T .Though this work is conducted in a flow regime different from that in the aforementioned investigations, it pioneers the concept of how different grain dynamics may affect the flow boundary condition by linking all the three velocity scales involved in χ and Ω.However, we recall that the μ w /f (r) relation is reported for the rapid flow regime, unlike how (1.1) and (1.2) describe a monotonic decay of μ w /f towards the creeping regime.Hence this work will focus on the two more recent models as an attempt to understand the degradation mechanisms for μ w .
Both (1.1) and (1.2) describe a monotonic decay of μ w /f with the strength of the associated cancelling mechanism, which suggests a correlation between the model variable χ and Ω. Equivalently, we may infer a hidden relation between the microscopic grain fluctuation and rotation dynamics.Hence this work aims to investigate experimentally the boundary stress of a confined simply sheared granular flow in view of the grain dynamics to broaden our understanding of the granular flow boundary friction coefficient.Section 2 describes the experiment set-up, while the image and data processing methods are presented in § 3, together with the measured bulk dynamics.Section 4 discusses the scaling relations between |ω| and T , and the measured data are used to consolidate a wall friction model.We conclude the work in § 5.

Experiment set-up and force measurement
Experiments were conducted in the annular shear cell in figure 1(a), identical to that used in Artoni et al. (2018).We used 2.1 kg of polyoxymethylene (POM) spheres of diameter D = 5.9 ± 0.1 mm and density ρ s = 1400 kg m −3 packed to height H = 15D.A top weight of mass M w = 0.22, 1.1 and 5.4 kg was added but allowed to move in the vertical direction to apply constant loading stress over its projection area S A .The confined spheres were driven by a bumpy bottom wall (see figure 1b) at different rates O = 5.9, 23.4 and 117.2 r.p.m., which combines with the three M w to achieve nine flow conditions, as summarized in table 1.
A six-axis force sensor (ATI Nano17) was mounted to a cut wall piece of the same curvature and material to the shear cell container (shown in figure 1c).The sensing module was installed to a 20 mm × 20 mm square window centred at different heights z = 15, 30, 50 and 70 mm from the base to measure the local boundary force.The force signals were averaged over time for three repeated measurements under each loading and driving condition to estimate the mean stress components.The normal stress σ rr in general follows the hydrostatic pressure profile σ rr (z) = M w g/S A + ρg(H − z), as reported in 969 A7-3 figure 4(a) of Artoni et al. (2018).The wall tangential stress components in the shearing and the vertical directions, σ rθ and σ rz , were used to calculate an effective wall friction coefficient by μ w = (σ rθ ) 2 + (σ rz ) 2 /|σ rr |.We follow Savage & Hutter (1991) to conduct a set of sliding table tests to measure the intrinsic grain-wall sliding friction coefficient f = 0.24 ± 0.01 (see the Appendix for more details).This value was used to scale μ w before we present the depth profile of μ W /f in figure 2(a).Degradation of μ w /f away from the moving base is observed when the loading is sufficient (M w ≥ 1.1) to engage sphere dynamics to the base motion.With the lowest M 1 , the spheres may remain loosely in contact with each other, and a nearly constant μ w /f was measured for z/D < 9 but starts to grow toward the top plate with μ w /f > 1, suggesting the top load effect.As we aim to study the wall friction weakening, the data of these three cases with the lowest M 1 will not be considered in the following analysis.

Image analysis and bulk dynamics
To understand the depth variation of μ w /f in view of (1.1) and (1.2), we obtain the needed velocity information from the individual grain dynamics.The sphere motion over the region boxed by the thick blue rectangle in figure 1(a) was recorded by a lateral high-speed camera (Phantom Miro 320S) at a rate of 24-3400 frames per second (FPS).A front LED lamp (EFFILux EFFISharp) was used to generate a bright reflection spot on each sphere surface to permit particle tracking, and two LED panels (FOTGA LED430) were placed aside to provide nearly uniform illumination.A snapshot of one recorded image is shown in figure 3(a).
The reflection spot on each sphere was nearly circular and can be located by the method of circular Hough transformation (Hough 1962) and found off-centred as the LED lamp was not aligned with the camera.We compared the bright spot positions to the manually located true centres in one image to evaluate a mean deviation distance and orientation, as r/D = 0.09 ± 6.1 × 10 −3 and θ = 1.15 ± 0.09 rad from the horizontal, with both standard deviations (STDs) lower than 7%.We exploited this nearly uniform deviation to offset the bright spot locations to present the sphere centres, as shown by the red plus signs in figure 3(a).Refer to Lin & Yang (2018) for more details.As the illumination and camera set-up remained unchanged, we applied the same offset strategy throughout the image processing routine.The nearest neighbour method was applied to track each sphere centre in consecutive images to evaluate instantaneous particle velocity vectors u pi (Yang & Huang 2016).
To measure grain rotation, a nearly circular small black marker of diameter (d = 1 ± 0.21 mm) was spray-painted on each sphere.We followed Lin & Yang (2018) to extend a search circle of radius D/2 from each sphere centre to locate the associated marker pixels using a greyscale threshold.We used the mean position of these identified marker pixels to present the marker centroid, as shown by the green points in figure 3(a).The same nearest neighbour method was applied to locate the same marker in the subsequent image.Under the assumption that the rotation axis remains unchanged over a short observation duration of 10 −2 s, we can fit a plane to the arc swept by five consecutive marker positions in the least squares sense.The normal vector of this fitted plane when forced to pass the instantaneous sphere centre found in the first instant specifies the desired rotation axis vector that can be normalized to unity magnitude.From this, we can extract the radius of gyration r pi for the marker-swept arc together with its tangential velocity V θ,pi , and the division of the two scalars determines the angular speed |ω pi | = V θ,pi /r pi .A further combination with the normalized rotation axis vector gives the desired angular velocity vector ω pi .More details regarding the measurement algorithms can be found in Lin & Yang (2018).
To estimate a bulk property from these instantaneous sphere dynamics, an averaging routine was implemented as described below.An averaging box that spans width L b = 6D .The instantaneous information of the spheres that fell in the averaging box over an observation period was collected to compute a temporal mean to present a quasi-steady bulk property at the mid-height of the average box.Together, we calculated the standard deviation of the bulk value over the observation period.We took a 5 s observation period if FPS ≤ 1000, and changed to a 1 s period if FPS > 1000.Next, we shifted the averaging box by 1D across the observation window in figure 1(a) to evaluate the depth profiles of bulk dynamic properties to be discussed in the following.First, the bulk streamwise velocity was evaluated as u (z) = i u p ,i (t)A i (t)/ i A i (t), using the instantaneous velocity component parallel to the base u p ,i (t) of the ith particle and its projection area in the averaging box A i (t).Here, A i (t)/ i A i (t) represents a weighting factor considering all the other particle information that fell within the averaging box at the same moment.The velocity perpendicular to the base, u z (z), nearly vanished across the depth for all flow conditions and hence is not considered.
where the velocity gradient at the averaging box centre z c was taken into account to offset the possible effect from the chosen averaging box height (Artoni & Richard 2015a).It is worth noting that the measurement errors of granular temperature could be discussed further in Xu, Reeves & Louge (2004).The scaled results, T /(OR) 2 , in figure 2(c) exhibit another sharp trend variation around the velocity transition height z/D = 10.Next, we examine the overall grain angular speed |ω| in figure 2(d) after scaling by O.As grain rotation is induced by unbalanced torque from its interaction with the neighbouring grains and the shear cell wall, it fluctuates more severely in time than the properties deduced from grain translation motion.Hence we present the expected value of ω pi collected over the same time interval as that used for the velocity.All the |ω|/O values exhibited an exponential decay across the shear zone, and a trend conversion into the creeping zone around the same z/D = 10, similar to that observed on T (z)/(OR) 2 .For both |ω|/O and T (z)/(OR) 2 , an increase of the top load caused a steeper decay in the shear zone and a stronger variation in the creeping zone, compared to that found for u (z).This implies a hidden link between T and |ω| that has never been reported, and suggests an association between the micro-events considered in the two μ w models.
Finally, we also examined the solid volume fraction profile φ(z, t) = i A i (t)/L b H b (not shown) to confirm a nearly consistent value at φ ≈ 0.7 across the creeping zone and most of the shear zone, while φ(z) drops near the moving base when z/D < 5 due to shear-induced dilatancy.The higher value of φ may be attributed to the ordering induced by the wall.However, the fact that the transition of u , T and |ω| occurs in the plateau of φ will support that the grain packing condition is not the primary cause for the phenomenon, but some other microscopic grain dynamics.

Scaling relations
First, we examine the D|ω|/2 as a function of T in figure 4(a), where the bulk variables are manipulated to possess identical dimension.A nearly linear correlation between the two variables is discovered in the creeping and fast shear zones.To provide further evidence, we also extracted the data from contact dynamics simulations of confined shear flows as presented in Artoni & Richard (2015a,b).The simulation results were processed in the same way as we analysed the experimental data and are shown by the grey solid diamonds in the plot of the shear zone in figure 4(a) to confirm a similar trend.Fitting to each data group on the log-log plot gives the best-fit slope, which diminishes slightly with the degree of streamwise fluctuation (≈ T ), from slopes of approximately 1. velocity at the contact point (i.e.|ω| D/2u 1) so that the grain-wall interaction is more like a pure sliding motion.Meanwhile, a small Ω suggests a large χ , which can also occur when the translation fluctuation is negligible as compared to the bulk motion (i.e.u / T 1).From both viewpoints, we may expect a μ w closer to the sliding friction coefficient between the grains and the wall.In contrast, when the strength of streamwise fluctuation increases (χ < 1), the intermittent grain motions are promoted to cause fluctuating interaction forces and unbalanced torque on each grain.The rotation index Ω rises accordingly to suggest how grain rotation diverts the sliding friction force at each contact to degrade bulk sliding.
In an attempt to link these micro-mechanisms to other flow variables, we recall that the streamwise flow fluctuation is often associated with the bulk shear strain rate | γ |, which is a pertinent variable in understanding fast flow behaviours.We evaluated | γ | from u (z) and examine how the two variables vary with | γ | in figure 5. We detect a clean power-law relation T ∼ | γ | 0.7 in the fast shear zone but the goodness of fit is lost in the creeping regime.The current power exponent in the shear zone is greater than the values, 0.5 to 0.625, reported from a simulated confined shear flow but losing the constant exponent in the creeping regime is a shared feature (Richard et al. 2020(Richard et al. , 2022)).Further, it is worth noting that our power exponent falls in the range of 0.5-1.0measured in the flowing layer of a quasi-two-dimensional narrow rotating drum flow (Orpe & Khakhar 2007) and below the exact 1.0 predicted for pure collisional flows in simple shear (Jenkins & Savage 1983;Campbell 1990).
The clean correlation between |ω| and T in figure 4(a) explains the similar variation trend for |ω| ∼ | γ | 0.84 in the fast shear zone and the scatterings in the creeping zone.To the best of our knowledge, how |ω| varies with γ at the lateral wall has never been reported.We may conclude that | γ | is no longer effective to characterize the particle level fluctuations in either their translation ( T ) or rotation (|ω|) in the creeping zone.Some may worry about the further influence from the top load or, equivalently, the confining normal stress, This inertial number measures the tendency of grain streamwise motion due to | γ | relative to its transverse settling under σ rr .Hence a flow at large I is in a well-fluidized fast state so that u can be large and the fluctuation-induced dynamic variations are comparably irrelevant in modifying bulk dynamics.This explains why figure 6 shows a general monotonic rise for χ but decay in Ω with increasing I.In particular, Ω saturates to a constant at large I to suggest a constant μ w /f from (1.2).Likewise, the discovery that χ rises rapidly with I also leads to a limiting constant μ w /f ≈ 1.These findings should explain why assigning a constant μ w works for fast flow prediction.
On the other hand, we notice a sudden rise of Ω for I 10 −2 -10 −1 , suggesting a variation in μ w /f .The particular inertial number I = 10 −2 has been marked to signal the transition from the quasi-static to the inertial flow regime in Midi (2004).In fact, Lu, Brodsky & Kavehpour (2007) also report a transitional regime in between the quasi-static and the inertial flow regimes over 10 −3 I 10 −1 (where I is converted here from the Savage number therein, which is exactly I 2 ).Hence we may expect the transition from the shear (fast inertial) flow regime to the creeping (transitional, towards the quasi-static) regime for the current flow condition.Finally, we would like to comment on the correlation between χ , Ω and I. Judging from figure 6(a), it is clear that the χ -I data do not collapse onto a universal trend, and that the Ω-I data show a pronounced scatter in the creeping regime in figure 6(b), both contrasting the universality of the Ω-χ relation in figure 4(b).Hence we would like to suggest that at least one of T and |ω| should be considered in a boundary condition model, in addition to the commonly used I.

Wall friction model
Finally, figure 7 examines the measured μ w /f with χ , and a monotonic decay with decreasing χ is observed, confirming the effect of how promoted fluctuations, relative to bulk u , can degrade the bulk friction coefficient from that of pure sliding.For the first time in the literature, we confirm in experiments that the μ w /f -χ model in (1.1) can capture the general trend of wall friction weakening with the fitted coefficients A = 0.1 and B = 2.05.The decay trend of μ w /f with decreasing χ is shown by the black solid line.However promising it is, one may be concerned that granular temperature is generally more difficult to measure than the shear strain rate, and may hope to express T by the fitted relation T ∼ | γ | 0.7 found in figure 5(a).However, such a nice correlation is not preserved in the creeping regime, where we detect the pronounced weakening of μ w /f .Hence the efforts seem critical and indispensable unless we may extract other robust correlations for T -| γ | in the creeping regime.
The other μ w /f -Ω model in (1.2) is proposed, with its two coefficients a and b representing the mean angular velocity components in the directions perpendicular and parallel to u .It would be meaningful to see how well this model performs while keeping the dependence on the just-confirmed parameter χ .Thus we further simplify (4.1) into Ω − χ −1 and rewrite (1.2) as A best fit to the experimental data gives ā = −0.1 and b = 0.39, and the new model shows equally good agreement, as portrayed by the grey dashed line.The negative ā indicates that ω ⊥ (∼ ā) points in the direction opposing u as grains rotate under the action of streamwise friction at the wall contact points.However, it is surprising to observe a positive b ∼ ω , which suggests that the individual grain-wall contact friction in the vertical direction did not cancel but showed a tendency towards the moving base even though the bulk transverse velocity is averaged to nearly zero.This phenomenon may be attributed to the grain-wall friction coefficient being much smaller than the grain internal friction coefficient f int = 0.394.This coefficient f int was estimated by the tangent of the static repose angle 21.5 • measured on a grain pile on the bed of the same experimental spheres.The grain-wall friction cannot compete with the inter-grain friction so that the grains adjacent to the container wall tended to rotate as sketched in the inset.
To support this speculation, we extracted the angular speed components of individual grains and used their temporal and spatial averages to estimate the bulk angular velocity components ω ⊥ and ω .The mean values ā = −0.2478and b = 0.3145 agree in sign with the above-fitted values, supporting the speculated mechanism further, and provide another pair of model coefficients.In addition, the standard deviations from the mean values a = ω ⊥ /|ω| and b = ω /|ω| measure the strength variation in the microscopic events and can be used to estimate a set of model coefficients.In figure 7, the green dotted line first portrays the predicted μ w /f trend using the mean (ā, b), which seems a bit off-trend with the measured data.However, if we take into account the deviation in the mean (ā, b), an upper bound of μ w /f is found taking the highest ā + a and the lowest b − b, as shown by the blue dotted line.Likewise, the orange dotted line gives the lower bound of μ w /f with the smallest ā − a and the largest b + b.It is interesting to note that the measured μ w /f -χ data are enveloped by these two bounding curves, supporting the proposed model and the current understanding that grain-level fluctuations -in either translation or rotation -are the key to degrading μ w from f .

Conclusions
This experimental work reports a non-constant lateral wall friction coefficient that weakens with the distance to the moving boundary when a loaded confined granular material was sheared at the base.To examine the data from the perspective of two existing models in (1.1) and (1.2), we measured the intrinsic sliding friction coefficient between the grains and the container wall, f , bulk shear velocity u , streamwise granular temperature T , and averaged grain rotation speed |ω|.The independence of the two models suggests a correlation between the model parameters χ = u / T and Ω = |ω| D/2u .In turn, we discover a hidden correlation between |ω| and T across the flow field, even though the two variables actually exhibit distinctive depth variations in the shear and creeping regimes (separated at z/D ≈ 10), just like u .This finding further helps to consolidate the two μ w /f models to propose (4.2) in terms of χ and physics-embedded model coefficients.Finally, we also examine χ and Ω with respect to the local inertial number I and observe a sharp rise of Ω when I 10 −2 in the creeping regime, suggesting that grain rotation or the equivalent translational fluctuations may bring non-negligible microscopic mechanisms to be considered in the creeping flow rheology.
To further test the robustness of the proposed model, the same confined shear flows but with different packing heights H are highly desired, as the thicknesses of the shear and creeping zones may change accordingly.Moreover, investigation in other Declaration of interests.The authors report no conflict of interest.

Appendix. Sliding table tests
To measure the intrinsic grain-wall sliding friction coefficient f , we conducted a sliding table test inspired by Savage & Hutter (1991).We designed an acrylic circular disk on which three experimental POM spheres were glued as the vertices of an equilateral triangle, as shown in figure 8(a).Special care was taken to ensure that the spheres were glued to the same height so that the disk can stand horizontally on these spheres.The disk was placed gently on an inclined plate made of poly(methyl methacrylate) (PMMA) in figure 8(b), identical to the material used for the cylinder wall of the shear cell.After we checked that the three spheres were in balanced contact with the plate, we slowly tilted the plate until the disk started to slide.At that moment, we recorded the inclined angle θ k by the electronic level meter so that the sliding friction coefficient can be estimated by f k = tan(θ k ).To reduce the possible effects from the set-ups, we alternated the orientation of the disk placement by positioning one of the glued spheres (labelled A, B, C) on the downstream side before the test.Once set in motion, the disk will slide as if the chosen sphere led the sliding down the plate, as illustrated in figure 8(b).We also alternated the plate placement by swapping the upstream and downstream sides.For each set-up, five repeated experiments were conducted and the obtained 30 f k values were averaged to give the mean value of 0.24 with a standard deviation of 0.01.

Figure 1 .
Figure 1.Experiment facilities.(a) Annular shear cell filled with POM spheres to a height H and confined by the top loading M w .The blue rectangle marks the observation window.(b) The base bumpy wall was rotated at rate O. (c) A force sensor was mounted on the lateral wall of the cell.

Figure 2 .
Figure 2. Depth profiles of bulk (a) scaled effective wall friction coefficient, (b) scaled streamwise velocity, (c) scaled streamwise granular temperature, and (d) scaled angular speed.A dashed horizontal line separates the shear zone (open symbol) and the creeping zone (filled symbol).

Figure 3 .
Figure 3. Portion of a raw image in greyscale.(a) The sphere (red plus) and marker (green asterisk) locations are marked.(b) The red and blue (moving 1D from the centre of the red box towards the bottom) boxes denote the averaging box used to extract the bulk properties.
Figure 2(b) shows the velocity depth profile u (z) after being scaled by the base rim tangential velocity (OR).Distinctive decay rates in the velocity depth profiles are observed.We follow Artoni et al. (2018) to define a fast-moving shear zone by the segment showing exponential decay to zero next to the moving base, and identify a much slower segment with a reduced decay rate in the creeping zone near the top plate.We fitted a line to each velocity segment and use their intersection to determine a transition height at z/D ≈ 9-11, as marked midway at z/D = 10.The velocity profile can be well described by u (z)/OR = u b + (u 0 − u b ) exp(−z/δ), as reported in Artoni et al. (2018).Here, u b and u 0 are the scaled velocities at the top and bottom, respectively, and δ is the decay length of the profile.The depth profile of the dominating streamwise granular temperature was evaluated next 969 A7-6 https://doi.org/10.1017/jfm.2023.

Figure 4 .
Figure 4. (a) Examination of the scaled angular speed D |ω|/2 with T in the creeping (inset) and fast shear zones.(b) Rotation index Ω versus χ.The contact dynamic simulation (CDM) data were adopted from Artoni & Richard (2015b), denoted by grey solid diamonds.

969Figure 7 .
Figure 7.Comparison of the experimental data of μ w /f -χ and the best fitted models in (1.1) and (4.2).Upper and lower bounds of μ w /f -χ are also portrayed, using the measured angular speeds as the model coefficients.

969Figure 8 .
Figure 8. Experimental set-up of the sliding table test.(a) The top and side views of the circular disk.Three spheres were inserted into the disk at positions A, B and C. (b) The sliding table: upstream and downstream positions are the initial placements of the disk.

Table 1 .
The details of the nine flow conditions considered in the experiments.