Characteristics of force chains in frictional interface during abrasive ﬂow machining based on discrete element method

Abrading is a very important sub-technology of the surface treatment technology with vast applications in the industry. This study aims at analyzing the inherent laws of friction systems during abrading. In particle ﬂow code modeling, the abrading process can be simpliﬁed to the movement of particles in a parallel-plate shear friction system. In this study, the PFC2D software is used to con-struct the particle ﬂow friction system with the set of parallel plates and the model parameters according to the abrading processing equipment and processing materials, control the simulation of a single variable, and compare the output data to estimate the impact of change of parameters on the force chain. The simulation results show that the shear dilatancy can be divided into three stages: plastic strain, macroscopic failure, and granular recombination stages. The distribution and load rates of the weak force chains depend on the load, velocity, friction coefﬁcient between granules, granular diameter, and number of granular layers. The number of granular layers and the load increase cause the direction of the force chain to be ori-ented with the vertical direction, and the force chains move toward the horizontal direction as the velocity increases. The increase in load does not cause the shear dilatancy stage to occur; the velocity, friction coefﬁcient between granules, and granular diameter increase cause the shear dilatancy to occur gradually.


Introduction
All granular systems, irrespective of whether conventional solid or a fluid matter, exhibit different properties because of the presence of particulate matter in the parallel plate of the granular flow frictional system. The shearing interface of the particulate matter with a granular media significantly affects many natural events and engineering applications, such as crustal movement, mud-rock flow, ore-rock crushing, grain processing, three-body friction, and abrasive flow machining [1][2][3][4][5][6][7][8][9][10]. The force chain networks produced by the granular systems under external load affect the movement and granular flow of the frictional systems, thereby affecting the quality of abrasive flow machining. This has attracted the attention of scholars globally.
The concept of particulate matter was proposed by French scientist Gennes, who is a 1991 Nobel laureate. In general, it refers to a large number of discrete particles forming a system, in which the particle size ranges between 1 lm and 2 lm [1]. According to the sparse degree of particle arrangement, a particle system can be divided into particle gas, particle fluid, and particle solid; the latter two belong to the dense particulate matter systems.
Bouchaud et al. [11] introduced the concept of force chains. The adjacent granules squeeze one another to form different force chain strengths under gravity or an external load. A force chain is the path on which the force passes along the granular contact network, while a contact network, in which some contact deformations are large and connected into a quasi-linear system, passes a significant force of gravity or a significant number of external loads to form strong force chains. The weak force chains are formed because of the weak contact deformation among granules and low external force. The strong force chains can bear load, and the direction of the weak force chains is perpendicular to the direction of the strong force chains to support the strong force chains and prevent breaking. However, changing the occlusion among granules enables the weak force chains to bear the load. The external load borne by the system is transmitted through the force chain network. The force chain network is formed by an external load acting on the particles that squeeze each other in the granular system. The force chains are extremely sensitive to changes in external load and geometric characteristics of the granular system. A slight change in the external load causes a large variation in the force chains of the same granular system [12]. Along with the formation of force chains, the geometric transition in force chain connectivity is also observed. The spatial distribution of force chains in the shear band exhibits two dominant length scales, i.e., force chains are separated by lateral support, on an average, by one particle diameter in one direction and two particle diameters in a second direction [13]. Gendelman et al. [14] presented a new formalism that employed the knowledge of the external forces and the orientations of contacts between particles (of any given size), to compute all interparticle forces. Using a geographical null model constrained by the particles' contact network, Bassett et al. [15] extracted chain-like structures, which were the reminiscent of force chains, and proposed three diagnostics to measure these chain-like structures.
It is difficult to provide an exact description of the force chains because of the complex disordered characteristics of the internal force chains in the granular system. Many researchers use different methods to analyze force chains under different conditions. Force chains were observed to play an important role in the microstructure of granular materials. In statics, the buckling of force chains is mainly located inside the shear band [16]. The International Fine Powder Research Institute (IFPRI) had funded an extensive research on dry powder and granular flows, especially dense flows at relatively high shear rates, in which the effect of force chains and jamming interactions were investigated considering the flow, stress and packing dynamics [17]. Besides, some researchers analyzed granular flows using simulation methods. Guo and Curtis [18] focused on the modeling of complex granular flows employing the discrete element method (DEM) approach, where the DEM models were applied to study the flow behavior of aspherical, flexible, or cohesive particles, including particle breakage. Azéma and Radjai [19] analyzed inertial granular flows and observed that contact anisotropy, force chain anisotropy, and friction mobilization were related to contact network and force transmission. Booth et al. [20] proposed that a continuum of friction angles existed between single-grain and bulk friction angles due to the grain-to-grain force chains. Marteau and Andrade [21] measured the force chains in opaque granular matter under shear using digital image correlation (DIC) and granular element method (GEM) to understand the micromechanical response of a complex granular assembly applied to macroscopic strains and stresses. Nicot et al. [22] highlighted the crucial interplay between force chains and adjoining clustered structures (grain clusters in three-dimensional (3D) conditions and grain loops in two-dimensional (2D) conditions). Further, to have a better understanding of force transition and stability mechanism of granular materials on mesoscale, a series of photoelastic experiments were conducted. The formation and evolution of force chains in the particle flow processes were also observed [23,24].
A traditional fluid theory cannot be applied to granules because of the complex nature of granular matter. However, a force chain reflects the internal evolution of the granules under external conditions on mesoscale. The evolution of granular matter in a granular flow frictional system can be explored through the study of the mesoscale force chains, and the evolution of granules can be described from a new aspect based on the simulation advantages of the DEM or finite element method (FEM) models. Tordesillas et al. [25,26] observed that the stability of force chains was highly correlated with the strength and volume properties of the granular material. The buckling of force chains resulted in noninterference and formation of shear bands in the granular material. Moreover, the mechanical response of dense granular materials indented through the rigid flat punch was closely related to the force chains and the surrounding particle clusters based on the subsequent studies. Meng et al. [27] studied the distribution of granular flow lubrication, evolution law of force chains, their dynamic state, and distribution probability of the force chain distribution characteristics by establishing the granular flow lubrication model on friction of parallel surfaces under sliding conditions based on the DEM model. Sun and Wang [28] observed that the direction of the strong force chains coincided with the direction of large external load through the discrete simulation of coplanar sphere centers of 12 000 of 2D equidistant granules, which were in static stacking. Wang et al. [29] established a parallel-plate shear cell to simulate the shearing of an infinite parallel plate. The change in the relevant parameters and the shear dilatancy can be divided into plastic strain, macroscopic failure, and granular recombination stages. Moreover, the force chain direction and the load behavior are described by the force chains ratio of the x component to the y component, load distribution curve, and force chains pattern. Before that, they studied the mechanical properties and the force chain variations of the granular flow in the friction pair gap to simulate the Taylor-Couette shear model using DEM [30]. Bai et al. [31] investigated the complex granular flow and mixing in a cylindrical bladed mixer (CBM) by the continuum approach of Eulerian FEM. The method is verified by comparing various features of the solid flows reproduced by FEM, such as the formation of heaps, recirculating flow around blades, and spatial distributions of velocities, along with those reproduced by DEM. Meng and Liu [32] analyzed the average velocity, velocity fluctuation, regional partition, and self-diffusion characteristics of the dense shear-granular flow in a parallel plate. Cai et al. [33] studied some interesting features of 2D granular shearing flow by molecular dynamic approach for a specific granular system. Wang et al. [34] simulated the Taylor-Couette shear model using DEM to study the mechanical properties and force chain changes in the granular flow in friction pair clearance. Chang et al. [35] designed and developed a friction equipment to study the frictional behavior between work piece surface and granular matter. Further, the velocity, attack angle, surface morphology, and particle size were considered in the experiments to reveal the effects of these input parameters on the force transmission behavior. Wang et al. [36] used a multi-scale method based on coupling the FEM with DEM to model the shear process of two parallel plates with solid particles between firstbody deformation and third-body rheology in a three-body interface to account for the interaction of them. Wang et al. [37] simulated the Taylor-Couette shear model using DEM to study mechanical properties and force chain changes in the granular flow in friction pair clearance. Meng et al. [38] analyzed the effects of upper plate friction coefficient on macroscopic granular flows, mesoscopic force distributions, macroscopic stress-strain responses, and evolutions of super-strong force chains using a discrete element analytical model that used the PFC2D software with an average solid fraction of 0.80.
Elkholy and Khonsari [39] used steel balls with a particle diameter of 3 mm as the third body, and shearing tests were performed using shear devices with different friction surface roughnesses. The shear force between the upper and lower plates and the displacement of the lower plate toward y axis were tested. The test showed that the displacement of the lower plate in the y axis was related to the surface roughness of the upper and lower plates, normal load applied to the upper plate, and rotation speed. They also compared the experimental results with simulation results and observed a good agreement between the two.
Zhang [40] established a parallel-plate shearing model using DEM to observe the movement of the particles between the upper and lower plates through the shearing motion of the upper and lower plates. The results of the study showed that the shear dilatancy caused by the threebody particles was the main cause of the frictional interface displacement.
Granular flow is a complex and enormous system. Many scholars have explored the characteristics of granular matter through a variety of experiments and methods from different aspects. However, many unique phenomena and mechanical properties of the granular matter are determined by the complex dynamic response of the force chain network, and very few studies on the factors that affect the force chain characteristics have been performed. In this study, a parallel plate granular flow system is constructed to simulate the force chains of granules based on the DEM model. The pattern of force chains is observed under the influence of load. The load and distribution rates of the weak force chains are analyzed based on the overall motion and shear dilatancy stages, respectively, as well as the direction of the strong-weak force chains to observe the change in the internal force chains characteristics of the granular system with its parameters.
Based on the aforementioned existing research work, this study further develops the parallel-plate granular flow model of a frictional system based on DEM, aiming to simulate the abrasive flow machining process, provide a method to analyze the force chains in granular flow, and clarify important relationships between the force chains and processing conditions. The remainder of this study is organized as follows. Section 2 briefly introduces the simulation model and related parameters. Section 3 is composed of four parts: part 1 analyzes the distribution and load rates of the weak force chains changing with parameters, such as load, speed, friction coefficient, and number of granular layers; changes in the force chain direction are studied in part 2; shape of force chains changing with load in different stages are analyzed in part 3; and finally, the shear dilatancy stage changing with various parameters are analyzed in part 4. Section 4 is the test verification of this study. Section 5 concludes this study.

Simulation model
Currently, the experimental study on the internal force chain of particulate matter indirectly explores the distribution of particle chains inside the force chain through the observation of the deformation caused by the particle material during extrusion using the photoelastic device. However, an experimental method is still not available for a force chain in the parallel-plate particle flow friction system, because the particulate matter is small, and the system is in the interior of the friction interface. Therefore, this study uses the DEM model to simulate and analyze the Characteristics of force chains in frictional interface during abrasive flow machining… 357 size and shape evolution of the force chains using the PFC software.

Physical model
In this study, a 2D Hertz-Mindlin contact model of the granular material Poiseuille shear flow is established using DEM for the ideal state of the granular flow, as shown in Fig. 1. The upper plate is the loading surface, and the loading surface of the model consists of minimal granules, which are exerted by the y axis direction pressure and can only move along the y axis. The driving surface is the lower plate and the driving surface moves uniformly along the x axis. Consequently, the periodic boundary condition is applied at the boundary of this model. The relevant parameters in the model are set according to a specific test, in which the granular material is silicon carbide. The loading surface material is brass, and the driving surface material is cast iron. The specific parameters of this model are presented in Table 1.
The simulation shows that there is a shear dilatancy in the particle system. Figure 2 shows that the location of the particles in the shear dilatancy process changes, the arrow representing the moving direction of the particles and the driving surface at this time. Particle A always moves along the x axis under the action of the driving surface. In Fig. 2a, A exerts a pressing thrust on C to move C to the upper right. In Fig. 2b, the particles are subjected to equilibrium forces along the y axis; C reaches the highest point and moves horizontally. In Fig. 2c, due to the continuous movement of A, the resistance of A to C decreases along the y axis. C moves in the lower right direction under the action of the load and the tangential force between A and C. As the shearing movement progresses, the process repeats itself. Therefore, only one shear dilatancy process is studied in this study.
In the process of shear dilatancy, due to the action of the driving surface on the particle system, relative movement between the particles occurs, which affects the loading surface and changes the displacement of the loading surface.  Poisson's ratio of the particle 0.15

Fig. 2 Demonstration of the shear dilatancy process
Furthermore, we can divide the process of shear dilatancy into three sub-phases: (a) plastic strain, (b) macroscopic failure, and (c) granular recombination stages (see Fig. 2).

Numerical simulation
DEM is a numerical method, proposed by Cundall and Strack [41] in the 1970s, to study the mechanical behavior and motion of discontinuous particulate matter. The basic idea of DEM is to treat the particulate matter as a set of rigid discrete elements. Initially, each granule is in equilibrium. When the external conditions change, some of the granules have a corresponding movement or displacement under gravity and external load, which contact with other granules resulting in a new state of force according to the contact model. More granules obtain this new movement or achieve displacement because of the transmission of the force. The force-displacement law correlates the forces and displacements between the granules that are in contact with each other or between the granules and the walls. The normal force component can be calculated according to where k n is the normal stiffness of the contact point, U i the normal contact displacement vector, and i the direction of the coordinate axis. The tangential contact force is calculated in increments, and the tangential contact force component obtained by superposition is given by Eqs. (2) and (3) where DU s i is the contact displacement increment, DF n i the tangential increment of the elastic contact force, k s the tangential stiffness of the contact point, and l the friction coefficient.
In this study, the Hertz-Mindlin contact model is a nonlinear contact model that is closely based on the Mindlin and Deresiewich and the Cundall models. The theory is defined by the shear modulus G and Poisson's ratio m of two spheres that are in contact with each other.
The normal stiffness of the contact is calculated by Eq. (4) and the tangential stiffness of the contact is calculated by Eq. (5) whereG is the equivalent elastic shear modulus,m the equivalent Poisson's ratio, andR the equivalent radius. The flow chart of the simulation is shown in Fig. 3. The requirement simulation model consists of three basic components: the composition of the granular matter; contact characteristics; and definition of the material properties, boundary conditions, and input of initial parameters.

Parameter description
The friction coefficient of the driving surface in the simulation process, that is, the ratio of the tangential and normal force l ¼ F s =F n is caused by the sliding and rolling frictions of the granules during the movement. A high friction coefficient denotes a high tangential force. The direction of the force tends to be along the axis, and more rolling friction can be observed on the drive surface. When the resultant force of a granule changes into its non- equilibrium state, the resultant force is the non-equilibrium force, and the maximum value of the non-equilibrium force of all the granules is the maximum non-equilibrium force, which is calculated by Eq. (6) where P F i ! m is the sum of the m forces of the ith granule, and the higher the value of F u is, the more severe the local motion will be in the particle system. The component ratio of the non-equilibrium force in the x and y axes can reflect the motion directionality of the granular system, which is calculated by Eq. (7) where F x and F y are the components of forces F on the x and y axes, respectively, and the larger the value of Q u is, the more biased the motion of the local granules is toward the y axis. The component ratio on the x and y axes of the force chains reflects the overall direction of the force chains, which is calculated by Eq. (8) where f x and f y are the components of the contact force f between the granules on the x and y axes, respectively, and the contact force is not zero. The larger the value of Q c is, the more biased the direction of the force chains is toward the y axis. Otherwise, it is biased toward the x axis.
The distribution and load rates of the strong-weak force chains between the granules reflect the distribution and loading capacity of the strong-weak force chains. The distinction between the strong and weak force chains is numerically based on the average force hf i, which is used as a critical value. When the value of the force chain is less than hf i, the force chains are weak. Otherwise, they are strong force chains.
Considering the weak force chains as an example, the distribution rate and load rate are calculated by Eqs. (9) and (10) where f i \hf i and i and j are the codes for the contact force between the granules. f i ð Þ 0 and f j À Á 0 indicate that the corresponding force chain appears only once. The distribution and load rates of the strong force chains are calculated by Eq. (11) 3 Results and discussion 3.1 Distribution and load rates of weak force chains changing with parameters The distribution and load rates of the weak force chains that change with load in two typical stages of the simulation are shown in Figs. 4a and b. The load and distribution rates of the weak force chains gradually decrease with an increase in load during shear dilatancy, as shown in Fig. 4a. The contact force chain network of granules in the shear dilatancy process, which makes the strong force chains increase and the supporting effect of the granular system stronger, is observed by a change in load. The distribution rate of the weak force chain is almost constant, and the load rate increases slightly with an increase in load in the whole motion stage of the granular system, as shown in Fig. 4b. That is, in the whole motion stage, the force contact network of the granular system is uniform, and the increase in load strengthens the weak force chains.
The distribution and load rates of the weak force chains that change with velocity in two typical stages are shown in Figs. 5a and b. The load and distribution rates of the weak force chains increase slightly with an increase in velocity, and the change in velocity has a little effect on the contact force chains between the granules in the shear dilatancy stage. In the whole motion stage of the granular system, the distribution rate of the weak force chains is uniform when the velocity increases, and the load rate decreases slightly. The force contact network of the granular system is uniform, and the increase in velocity reduces the strength of the weak force chains.
The distribution and load rates of the weak force chains vary with the friction coefficients, as shown in Fig. 6. When the friction coefficients are 0.21, 0.23, 0.35, and 0.37, the granular system is in the whole motion stage. When the friction coefficient is small, the distribution and load rates of the weak force chains decrease with an increase in the friction coefficient. The contact among the force chains is even, which decreases the distribution rate of the weak force chains and increases that of the strong force chains because of the enhanced interaction between the granules. The distribution and load rates of weak force chains increase as the friction coefficient increases. When the friction coefficient is high, the squeezing between granules increases the contact between the weak force chains because of the resistance of the force chains. When the friction coefficient changes from 0.25 to 0.35, the change trend of the weak force chains is different because of shear dilatancy. The load rate of the weak force chains increases slightly. However, the distribution rate is basically unchanged. The increase in the tangential force of the driving surface increases the weak force chains slightly, while the distribution of the force chains of the granular system remains the same.
The distribution and load rates of the weak force chains vary with the number of granular layers, as shown in Fig. 7. First, the distribution and load rates of the weak force chains in the granular system increase and then decrease gradually. When the number of granular layers is three, the force chains form a whole, which connects like a network under the normal and tangential forces of the   driving surface because of the less number of granular layers. The distribution of strong force chains is wide and the local independent strong force chains are few, which make the distribution and load rates of weak force chains low. The amount of elastic compression in the granular system as well as the force chain contact among the granules increases when the number of layers increases. The total amount of elastic compression between granules continues to increase with the number of granular layers, so that the original micro-deformation contact, which subsequently forms into strong force chains, results in the reduction of distribution and load rates of the weak force chains, to strengthen them gradually.
The distribution and load rates of the weak force chains vary with the granular diameter, as shown in Fig. 8. The distribution and load rates of the weak force chains are small. The distribution rate does not change with the increase in granular diameter. The granular system is in the whole motion stage when the granular diameter is small, resulting in low distribution and load rates. Moreover, an increase in granular diameter has little effect on the distribution and load rates during the shear dilatancy process.

Force chain direction analysis
The whole direction of the force chain gradually changes from the x axis to the y axis with the increase in load, as shown in Fig. 9. In the shear dilatancy stage, the direction of the force chain changes more than the whole motion stage. The anisotropy of the granular system is changed so that the contact force on the x axis is reduced, and to avoid the shear dilatancy between the granular layers that occur because of the change in the whole direction of the force chain inside the granular system. The change in the whole force chain direction is related to the change in the friction coefficient on the driving surface. The friction coefficient decreases step by step and the tangential force along the x axis decreases gradually, so that the force chain from the drive surface to the inside of the granular system weakens along the x axis as the load changes, thereby affecting the whole force chain direction.
The force chain direction of the granular system in the shear dilatancy and whole motion stages changes with velocity in the simulation, as shown in Figs. 10a and b. The overall force chain direction gradually shifts from the y axis to the x axis as the velocity increases. Moreover, the force chain direction in the shear dilatancy stage changes more than in the whole motion stage. The direction of the weak chain is in an overall fluctuation state, in which the change in direction with change in speeds is relatively small. The change in the overall force chain direction is related to the change in the friction coefficient on the driving surface. That is, the friction coefficient increases gradually with speed changes. Thus, the tangential force along the x axis increases gradually, and the force chains from the driving surface to the inside of the granular system increases along the x axis, which affects the overall force chain direction.
The force chain direction changes with the number of granular layers, as shown in Fig. 11. When the number of layers is three, the direction of the weak force chain is more inclined toward the y axis than the strong force chain. However, with an increase in granular layers, the strong force chains are more biased toward the y axis than the weak force chains. The component ratio of force chains is greater than zero, that is, the overall direction begins to be biased toward the y axis. The distribution of the overall tilt grid force chain of the granular system is wide and the force chains are holistic, which are more along the x axis because of the effect of the driving force. The integrity of the granules is weakened with the increase in granular  layers, and the tangential force cannot make the force chains along the x axis between the upper granules. Thus, the strong force chains along the y axis can continue to exist under the load.
The force chain direction changes with friction coefficients, as shown in Fig. 12. In the whole motion stage of the granular system with the friction coefficient of 0.21 and 0.23, the strong and weak force chains are obviously separated, and the weak force chains are biased toward the y axis and the strong force chains are biased toward the x axis. In the whole motion stage of the granular system with the friction coefficient of 0.35 and 0.37, the weak force chains evolve toward the x axis with the increase of the friction coefficient. It is more biased toward the x axis than the strong force chains, and the strong force chains evolve toward the y axis. Due to the increase of the internal friction coefficient of the granular system, the structure of   force chains in the granular system changes significantly, the interlocking between the granular systems increases and the strong force chains along the y axis can exist with the interlocking.

Component ratio of force chains
The direction of the strong and weak force chains of the granular system, based on the simulation, is evident at the intersection A when the load is 2.56 9 10 -3 N. The whole motion stage of the granular system ends after the intersection, and shear dilatancy occurs. The force chain direction is the same after the intersection of the strong and weak force chains component ratios and the contact force drives the granules motion, as shown in Fig. 13. The intersection of the component ratio is the critical point at which the force chain mode changes. Therefore, the component ratio of the force chains of the granular system is calculated before and after the intersection point, as shown in Fig. 14. The component ratio of the weak and strong chains increases gradually with an increase in load. When the load is P\2:56 Â 10 À3 N, the component ratio of the weak and strong force chains is lower than zero. Hence, the force chain direction is biased toward the x axis. The effect of the strong and weak force chains on the granular system is similar at the beginning, driving the granular system to move along the x axis, whereas the partial weak force chains support the strong force chains on the y axis. At this stage, the growth ratio of the weak force chains component is high, the growth ratio of the strong chains component is low, and the direction of the weak force chains and the effect of the force chains changes gradually. When the load is P [ 2:56 Â 10 À3 N, the component ratio of the weak force chains is greater than zero, and the component ratio of the strong force chains is lower than zero. The effect of the force chains is separated, and the weak force chains begin to support the strong force chains along the y axis. The strong force chain still drives the granular system. The role of the strong and weak force chains is similar to the static granular system, that is, the weak force chains support the strong force chains that bear the load. A deviation in force direction is observed. Therefore, when the load is P\2:56 Â 10 À3 N, the displacement of the loading surface is almost constant. The bottom granules promote the movement of the upper granules, resulting in the occurrence of shear dilatancy delay. Thus, no relative movement between the granules of the whole system can be observed, which is similar to that in a static granules system. The component ratio of the strong and weak force chain intersection increases gradually with an increase in load, and the shear dilatancy phenomenon does not appear at this time. From the contact mode of the granular system, the location of granules, whether they are in the shear dilatancy or whole motion stages, can be determined by the result of the interaction   between the tangential forces along the x axis and the support forces along the y axis.
The component ratio of the strong and weak force chains, at the time of the intersection and before, changes with the velocity, as shown in Fig. 15. The component ratio of the weak and strong force chains is slightly lower with a change in velocity. In comparison, the weak force chains change greatly. The component ratio of the weak force chains is higher than zero and that of the strong force chains is lower than zero.
At this time, the effect of the force chains is different, and the weak force chains begin to support the strong force chains along the y axis. The strong force chains still drive the granular system. The appearance of the whole motion stage is clear at this moment. After v = 0.4 m/s, the component ratios of the weak and strong force chains are lower than zero. The force chain direction is on the x axis, and the effect of the strong and weak force chains on the granular system is similar at the initial time. Shear dilatancy occurs continuously.
The evolution of the strong and weak force chain directions during the simulation process increases with time in five layers, as shown in Fig. 16. The force chain direction, before and after the intersection of the strong and weak force chain curves, is different. Before the intersection, the weak force chains are biased toward the x axis, and the strong chains are biased toward the y axis. After the intersection, the relationship between the two force chains is opposite. The effect of a tangential force on the upper granular layer is weakened because the granular system is in the shear dilatancy stage before the intersection; thus, the weak force chains increase along the x axis and the strong force chains along the y axis. In the macroscopic failure and granular recombination stages, the strong and weak force chains separated; the granular system was loose; the response to the external force was not obvious. In the macroscopic failure stage, as the entire granular system is in an elastic potential energy release phase, the component ratio of force chains rises; thus, the component ratio of the force chains changes a lot when t = 2.0 ms, and there is a clear separation between the component ratio of strong force chains and weak force chains. The weak force chains are biased toward the y axis to support the strong force chains, while the strong force chains are biased toward the x axis to drive the movement of the underlying granules.
The component ratio of the strong and weak force chains, before and after the intersection, changes with the number of granular layers, as shown in Fig. 17. Before the intersection, the direction deviation of the strong and weak force chains is gradually reduced to be almost unchanged, whereas after the intersection, the direction deviation of the strong and weak force chains is gradually reduced until they are unified.

Shape of force chains changing with load in different stages
The morphological evolution and unique characteristics of the force chains are analyzed at different shear dilatancy stages. The black lines represent the force chains. The larger the width of the line is, the greater is the value of the force chains. However, only the relative size of the force chain is applicable at the same time.
The force chains are widely distributed in the granular system from the driving surface to the loading surface, and the force chains are continuous in the initial stage. As shown in Fig. 18, on the plastic strain stage, the shape of the force chains is similar to that of the grid; thus the force chains are defined as the ''gradient grid'' force chains. The diagonal force chains, such as the force chains in Region C, appear in the ''gradient grid'' force chains in Region A with a change in load. The initial force chain directions are shifted to the y axis because of the increase in load, such that the diagonal force chains appear in the upper layer granules in Region C to withstand the force in the y axis. In Region B, a similar ''Y''-shaped force chain significantly increases with load, such that a part of the underlying granules have the same trajectory movement and no extrusion contact force between the granules occur, which indicates that the increase in load makes the contact between the underlying granules and the drive surface uniform The weak force chains along the y axis gradually disappear with time. As shown in Fig. 19, in the macroscopic failure stage, the force chain distribution between the driving and loading surfaces is scattered and discontinuous. The presence of the strong force chains between the underlying granules is not only related to the interaction between the driving and loading surfaces, but also to the inertial collision between the granules. The shape of the force chains is similar to a tadpole. Hence, the force chains are defined as ''tadpole'' force chains. In Region E, the driving surface drives the bottom granules to the right, which is the direction point; whereas the force chains in the upper layer of the granules move toward the lower left, which are caused by the load extrusion from the loading surface and the counter-clockwise rotation of the bottom granules. The force chain in Region F moves toward the upper right. At this time, the bottom granules of Region E are transmitted through the collision of the tadpole-like force chains; thus, the rotation rate is lower. The friction caused by the rotation between the upper granules is not obvious, and the tangential force of the driving surface is transferred to other working surfaces. When the load is P = 4.48 9 10 -3 N, the form of the force chains is distributed along the y axis in Region G that corresponds to Region E, and the ''Y''-shaped force chains increase gradually with an increase in load. As the load increases and the friction coefficient decreases, the extrusion between the upper and lower parallel plates increases, and the ''Y''-shaped force chains between the granules increase to support the pressure of the granular system and withstand a stronger load. At the end of the shear dilatancy process, the shape of the force chains is ''pinnate''. As shown in Fig. 20, the granules on the left side of Region H pass through the force chains of Region H, and the granular rate is lowered so that the extrusion force is generated between the upper and bottom granules to form the force chains that are distributed along the y axis. The length of the force chains in Region H decreases gradually with an increase in load. The force chain distribution in Region I increases with an increase in load. The granular system is closely arranged to one another as the load increases to limit the free movement of individual granules in the system. The force chains caused by the velocity collision are reduced, and the load makes the contact between the bottom particle and the driving surface uniform, thereby increasing the contact force and force chain distribution area.

Shear dilatancy stage changing with parameters
There are two indicators of change in the shear dilatancy stage: the number of the shear dilatancy stage, which is the number of stage that occurred within 0.02 s; and the duration of the shear dilatancy stage, which is the duration of a single stage. The force chain is continuously distributed from the driving surface to the loading surface during the shear dilatancy stage, especially in the initial stage, and the force generated by the driving surface can operate on the loading surface well. Therefore, in certain period of time (0.02 s), the longer the shear dilatancy lasts, the higher the efficiency of the abrading process is, which is more conducive to the control of the processing Fig. 19 Shape of force chains changing with load in the macroscopic failure stage Characteristics of force chains in frictional interface during abrasive flow machining… 367 conditions. As shown in Fig. 21, the parameter changes have little effect on the ratio of plastic strain time to total time. The number of granular layers parameter is an exception to this. When it changes, the ratio of plastic strain time to total time is decreased, which deteriorates the quality of abrading. The number and duration of the shear dilatancy stage that change with load in the simulation are shown in Figs. 22a and b. First, the number of the shear dilatancy stage increases and then decreases until the shear dilatancy stage occurs randomly or rarely with an increase in load. Thus, the shear dilatancy stage will not occur continuously and the proportion of the whole motion stage will increase. This will improve the quality of abrading.
The number and duration of the shear dilatancy stage that change with velocity are shown in Figs. 23a and b. First, the number of shear dilatancy stage increases and then gradually stabilizes with an increase in velocity. Thus, the proportion of the whole motion stage during the granular system decreases gradually, and the shear dilatancy stage occurs continuously until the number of the shear dilatancy stage reaches saturation. This will improve the quality of abrading.
The number and duration of the shear dilatancy stage that change with the friction coefficient between granules in the simulation are shown in Figs. 24a and b. The number of the shear dilatancy stage does not increase rapidly; on the contrary, it gradually decreases until the shear dilatancy disappears; while the duration of the corresponding single shear dilatancy stage increases with an increase in the friction coefficient. Thus, the whole motion stage of the granular system disappears; the shear dilatancy stage occurs; the number of shear dilatancy gradually decreases until it disappears; and the duration of shear dilatancy increases gradually. This will improve the quality of abrading. Figures 25a and b show the number and duration of the shear dilatancy stage changing with the granular layers in the simulation. The number of the shear dilatancy stage increases gradually with an increase in the number of granular layers, while the duration of the single shear dilatancy stage decreases, and the ratio of plastic strain time to total time also decreases. This will deteriorate the quality of abrading. Figures 26a and b show the number and duration of the shear dilatancy stage changing with the granular diameter in the simulation. The number of the shear dilatancy stage starts from scratch and then decreases gradually; and the duration of a single shear dilatancy stage gradually increases with an increase in the granular diameter, which means that the whole motion stage of the granular system disappears. This will improve the quality of abrading.

Experimental test
This study verifies the simulation result by setting a test and observing the changes in the gap value using a rotary polishing testing machine to achieve infinitely long shearing motion. We fixed the rectangular samples of a specific size on the machine, imported the particles, rotated the abrasive disc with sandpaper to abrade the rectangular samples, and used an eddy current sensor to collect the gap value data. About the processing condition parameters, the load value is set to four values of 7 500 Pa, 10 000 Pa, 15 000 Pa, and 20 000 Pa; for the speed values of 60 r/min, 90 r/min, 120 r/min and 150 r/min; and three particle sizes of the same material particles: 1# (0.2 -0.3 mm), 2# (0.4 -0.5 mm), and 3# (0.6 -0.8 mm).

Experimental scheme
According to the particle size, load and speed, the test plan is presented in Table 2. Figure 27 shows that the load has a little effect on the gap value. This may be because the interaction between the strong and weak force chains determines that the gap value will not change much when the speed is within a normal range (no excessive centrifugal force) and the particle diameter is the same. The fluctuation and decrease in the gap value of 2# and 3# particles may be due to the decrease in the load rate and the distribution rate of the weak force chains. At the same time, the increase in the load causes the

Gap values changing with speed
The relationship between the gap value and particle size is shown in Fig. 28. From Fig. 28, it can be seen that regardless of the processing conditions of the experimental scheme, the particle size has a significant effect on the gap value, and the gap value increases as the particle size becomes larger. At the same time, there are obvious differences in the gap value of the different particle sizes under the same experimental conditions, which means that the larger the particle size, the more significant the effect of load on the gap value. The increase of gap value may be related to the increase in the distribution rate of the weak force chains-the increase of the gap value leads to an increase in strong chains, which requires the support of more weak force chains.

Gap values changing with particle size
During the test preparation, it was observed that when the particle size was large, and the rotation speed was high, the silicon carbide abrasive paper used to protect the sample on the abrasive disk wears sharply and breaks in a short duration, resulting in a large error. Therefore, in the experiment of the interfacial gap characteristics at different speeds, only the 1# particles were used to perform at loads of 7 500 Pa and 20 000 Pa.  The relationship between the gap value and speed is shown in Fig. 29. It can be seen that when the speed is 150 r/min, the gap value is obviously reduced. It is observed that the turning point occurred due to the excessive centrifugal force generated by an increase in speed. At the same time, the difference in the turning point of 7 500 Pa (90 r/min) and 20 000 Pa (120 r/min) indicates that the load is the influencing factor of the third body particle importing. The distribution and load rate of the weak force chains, which do not change much, may cause a small change in gap value before 150 r/min. At the same time, a decrease in the gap value caused by excessive centrifugal force also reflects the phenomenon that the direction of the force chain is biased toward the x axis when the speed is increased, indicating that the speed during the processing should not be too big.

Conclusions
It is important to explore the evolution of the internal force chains in the friction system between the parallel plates, because it has a great significance in the research on how to improve the surface quality at abrasive flow machining. In this study, a 2D Hertz-Mindlin contact model of the granular material Poiseuille shear flow is established using DEM for an ideal state of the parallel-plate granular flow. The direction of the strong and weak force chains is investigated using external and internal parameters as well as the distribution and load characteristics of the weak force chains. The external parameters include the load and the drive surface speeds, and the internal parameters include the friction coefficient between the granules, number of granular layers, and granular diameter. The load and distribution rates of the weak force chains are affected differently by those parameters. The intersection of the   strong and weak force chain component ratios can be used as a critical point in the whole motion stage. The greater the value of the intersection point is, the more likely the whole motion will occur. An increase in load and number of layers cause the force chain's direction to be biased toward the y axis, and an increase in velocity causes the force chain direction to be biased toward the x axis after v = 0.4 m/s. The shapes of the force chain have different changes under the influence of different loads. The force chain distribution is concentrated in the initial stage of shear dilatancy, and the abrading process has a high efficiency. It is easy to control the abrasive machining precision. The simulation results show that an increase in load and granular layer will deteriorate the quality of abrading, and an increase in velocity, friction coefficient between granules, and granular diameter will improve the quality of abrading. The experimental results indicate that the gap value does not change much with the processing conditions (except the particle diameter, which directly affects the gap value), and the processing speed should not be too high. These are also in agreement with the simulation results up to some extent.